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ABSTRACT 


Cross body thrusters permit a body of revolution Autonomous Underwater Vehicle to 
retain the energy efficiency of forward travel while increasing the ability to maneuver in 
confined areas such as harbors and piers. This maneuverability also permits more 
deliberate underwater surveys using a fixed, mounted forward and downward looking 
sonar. This work develops the necessary hydrodynamic coefficients, using methods 
applied to earlier vehicles, to develop a valid computer simulation model. Additionally, 
this work develops a polynomial regression translating thruster input in RPM to an 
applied force output, which is incorporated into the vehicle model. This model is then 
employed to examine the response and control, specifically at low speed, of a body-of- 
revolution Autonomous Underwater Vehicle equipped with off-axis cross-body thrusters. 
These results are then utilized to develop a series of PID controllers for use onboard the 
REMUS Autonomous Underwater Vehicle. 
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I. INTRODUCTION 


A. AUTONOMOUS UNDERWATER VEHICLES — A HISTORICAL 
PERSPECTIVE 

The concept of an underwater vehicle is not new; Leonardo da Vinci envisioned 
submersible troop transports, assault craft, and diving rigs in the early sixteenth century. 



Figure 1. Da Vinci’s sketch of a manned submersible, upper-left, [1] and diving rig, 

right, [2] circa 1500. 


Ultimately, the idea remained on paper only, as full realization extrapolated from 
practical test results was deemed too dangerous by Da Vinci [3]. Nearly three centuries 
later, the Bushnell brothers built the first American submarine, the Turtle, in 1775. The 
tiny, single-seat, egg-shaped craft attempted the first recorded, and ultimately 
unsuccessful, submarine attack in New York Harbor in 1776, against the HMS Eagle [4]. 
Military applications of underwater vehicles have steadily marched forward since. 
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Parallel to the development of manned undersea vehicles, came the creation of undersea 
weapons systems, principally the advent of the torpedo—the first true autonomous 
underwater vehicle [5]. 

The Naval Warfare strategy of the United States constantly evolves and adapts to 
a varied array of threats in the dynamic and globalized modern world. Technological 
advances in recent years have seen a tremendous increase in the use of autonomous or 
remotely operated vehicles. These vehicles provide combat commanders a real-time 
intelligence and communications capability previously unheard of on the battlefield. As 
the focus of naval operations has shifted from the open ocean (“blue water”) to the littoral 
(“green water”) zones, the need for a capable and advanced autonomous underwater 
vehicle (AUV) has steadily increased. 

The AUV provides commanders a vital near-shore intelligence and 
communications asset. For example, a mission profile may require the ingress of special 
operations personnel onto a beach via submarine. The AUV would enable the mission 
commander to receive real-time data from the proposed line of advance, without 
endangering any personnel. Furthermore, that same AUV could be used as a 
communications relay from the special operations team back to the submarine. 

In the civilian sector, AUVs and their tethered brethren, remotely operated 
vehicles (ROV), have been used extensively in salvage operations. Woods Hole 
Oceanographic Institution’s ROV Argo, for instance, was a primary tool used in 
discovering the wreck of the Titanic [6]. The tethered submersible enabled 24-hour 
reconnaissance of the wreck site at depths of almost 2.5 miles while the research teams 
remained above the surface. 

More recently, the Deepwater Horizon (British Petroleum) Oil Spill in the Gulf of 
Mexico, employed teams of ROVs in several efforts to stem the flow from the leaking 
pipe. These ROVs featured fully maneuverable actuator arms and were operating in 
depths of nearly 1 mile [7]. These depths are well beyond the limit of human physiology, 
as the depth limit for a surface supplied diver is 285 feet [8]. 
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The AUV provides the next step in the advancement of underwater robotic 
technology. While the removal of the human element entirely from the vehicle control 
interface requires complex and detailed control algorithms, the benefits far outstrip the 
risks involved in unmanned undersea vehicle operations: 

• Un-tethered vehicles enable reconnaissance through a larger range, 
providing a cost-effective method of narrowing a search area for future 
study 

• The AUV is programmed with its own logic algorithms, enabling 
instantaneous on-board decision making based on current sensory 
information 

• No external interface is required - that is, all power and maneuvering 
controls are internal to the vehicle 

• The AUV can operate in depths and temperatures human divers cannot 

• Generally, AUVs are cheaper than ROVs, as they do not require a tethered 
support platform with the requisite human interface for mission 
accomplishment 

B. MOTIVATION 

Of specific interest to this thesis, is the use of AUVs in mine countermeasure and 
clearance operations. AUVs provide a safe and cost-effective means of determining the 
size of a minefield and the specific locations of individual mines. Future AUVs may 
even incorporate onboard mine countermeasures to defeat or disable encountered mines 
without endangering Naval personnel; providing a similar capability as the remotely 
operated ordinance disposal robots provide to U.S. military and law enforcement 
departments. To facilitate this, precise control in all six degrees-of-freedom during low- 
speed operations is required. This work will study the control and guidance of a body-of- 
revolution (“torpedo-like”) AUV fitted with cross-body thrusters. These thrusters 
provide greater control and maneuverability at slow speeds than a standard single 
propulsor counterpart. Control surfaces require a minimum speed of advance in order to 
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provide effective lift to the body; consequently, this speed may be faster than the optimal 
resolution capability of the onboard sensors, specifically the sonar. Cross-body thrusters 
provide near-instantaneous control in all six degrees of freedom, enabling full 360-degree 
analysis of an “interesting” object as well as more “cost-effective” path planning, 
waypoint navigation and hovering methods. 

C. LITERATURE REVIEW 

Cross-body thruster control of an AUV is not a novel idea. Several AUVs, past 
and present, use cross-body thrusters including, among others: NPS Phoenix (c. 1996— 
1999), NPS Aries (c. 2002), REMUS (c. 2001-present) and Odyssey IV (2009-present). 



Figure 2. AUVs past and present - NPS Phoenix, top left [9], NPS ARIES, top center 
[10], MIT Odyssey IV, top right [11], NPS REMUS, bottom 

The initial design of the REMUS underwater vehicle did not feature cross-body 
thrusters and was instead controlled by a main propeller astern and a cruciform series of 
control fins. Prestero [12] studied this vehicle and proposed both non-linear and linear 
models for motions in all six degrees-of-freedom based on rigid-body Newtonian 
physics; that is: 


y^F = m-a 

(1) 

YjM=I-™ 

(2) 


Model linearization is standard and provides a simple way to de-couple the otherwise 
highly non-linear equations of motion. The key tenet to this solution is its reliance on 
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small motions and small changes in the system dynamics; thereby eliminating several 
non-linear cross terms, keeping damping terms quadratic, and diagonalizing the mass and 
inertia tensors. Furthermore, restricting the model to small perturbations allows the 
matrix of Coriolis terms to be neglected. 


Prestero did voluminous work in determining the exact parameters and 
coefficients of the original REMUS vehicle in operation; necessary for complete 
modeling of the vehicle in all six degrees-of-freedom. The values were detennined 
analytically using several methods and then compared to data detennined from 
experiments using the actual REMUS vehicle in the tow-ta nk . While several of the 
coefficients, including fin lift coefficients and forces, will remain close to the original 
values; the current evolution of the REMUS vehicle will have several changes to the 
added mass coefficients; as several of the principle dimensions are significantly larger. 
This thesis will utilize several of the same methods in determining the hydrodynamic 
coefficients and parameters pertaining to the enhanced version of the REMUS vehicle, 
using the original parameters determined in [12] as a baseline for comparison. 


Fodrea [13] focused on utilizing forward-looking sonar in path planning and 
object avoidance maneuvers, using a REMUS vehicle of the same series as Prestero [12]. 
This problem has been under research since the advent of the autonomous underwater 
vehicle [13] as path planning and object avoidance are critical to mission success and 
vehicle longevity. The work proposed a multi-variable sliding mode controller (SMC), 
particularly useful in linear systems, as the governing control method for the line-of-sight 
guidance algorithms. The sliding mode controller provides a robust control system, 
advantageous to autonomous vehicles in the underwater environment. Fodrea reduced 
the system to two-dimensions (x-y plane only). The line-of-sight guidance algorithm 
functions by attempting to drive the heading error to zero: 


where 


Vws^Vu-aciM-vit) 


U track (0 = tan 




( 3 ) 

( 4 ) 
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The work resulted in a very responsive controller, allowing the REMUS vehicle to regain 
desired track immediately after avoiding an object [13]. Crucial to this was the design of 
the weighting functions for bearing and range. The forward-looking sonar provides a 
very specific cone of view, and its accuracy is highly subject to the speed of the vehicle. 
The weighting functions applied large (“heavy”) values to objects close to the vehicle, 
decreasing in value as range and relative bearing from the centerline increased. 

The thruster control problem in a shallow water environment has been 
investigated by Reidel [14] using the NPS Phoenix AUV. Riedel examined the effects of 
environmental disturbances (waves) on a low-speed vehicle, and developed a Disturbance 
Compensation Controller (DCC) to assist the NPS Phoenix AUV in maintaining station, 
or loitering, in a particular area of interest. The DCC employs an Extended Kalman 
Filter (EKF) to continuously update the estimated system states based on corrupted 
(“noisy”) measurements. This controller resulted in a new generalized approach to 
modeling underwater vehicles exposed to shallow water waves and currents, 
implementing direct fluid measurements to assimilate short-term wave-magnitude and 
direction infonnation into a general prediction of sea disturbances. 

Loitering control was examined further by Cooney [15], studying the MIT 
Odyssey IV AUV, which operated using azimuthing thrusters. The work presented a 
simple model for the unsteady thrust imparted by the azimuthing thrusters and developed 
a non-linear Model Predictive Controller to provide superior response when utilizing the 
thrusters in comparison to a standard Proportional-Integral-Derivative (PID) controller. 

D. SCOPE OF THIS WORK 

This work builds upon the foundation of REMUS knowledge previously 
presented, providing a robust vehicle simulator responsive in all 6 degrees-of-freedom. 
While the REMUS vehicle is highly responsive when using the cross-body thrusters, a 
valid simulation thereof remains unknown. This thesis develops a non-linear SIMULINK 
model in 6 degrees-of-freedom, modeling lift, drag, added mass, control forces and 
moments as an array of hydrodynamic coefficients. The non-linear model is then 
linearized to conduct analysis on lateral and vertical thruster effects on vehicle motion to 
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validate simulation results against actual REMUS vehicle in-water response and to 
develop a PID controller to minimize vehicle motions in the horizontal plane. The 
desired end-state is the implementation of this robust vehicle model and controller into 
the vehicle itself for further study. 

Chapter II provides general background information on the history and 
development of the REMUS AUV. This will include comparative tables of physical and 
functional characteristics between the initial vehicle studied in [12] and [13] and the 
current vehicle in possession at NPS. 

Chapter III presents the necessary vehicle and environmental assumptions, 
coordinate frame relations and equations of motion to model the REMUS vehicle in six 
degrees-of-freedom. 

Chapter IV derives the necessary added mass matrix, hydrodynamic coefficients, 
and thruster force and torque relations utilizing several methods presented in [12]. 

Chapter V combines the equations of motion determined in Chapter III with the 
complete hydrodynamic coefficients presented in Chapter IV to obtain the full matrix 
representation of REMUS vehicle motion across all six degrees-of-freedom. 

Chapter VI provides focused study on implementing control laws to eliminate 
vehicle instabilities as well as examining the full controllability of the vehicle. This 
chapter also provides an analysis of thruster effects on the linearized vehicle models. 

Chapter VII discusses options to reduce error and provides possible avenues for 
continuing research. 


7 



THIS PAGE INTENTIONALLY LEFT BLANK 


8 



II. THE REMUS AUTONOMOUS UNDERWATER VEHICLE 


A. HISTORY 

Developed in 1994 by von Alt and associates of the Woods Hole Oceanographic 
Institution, the demote Environmental Measuring Units, or REMUS, autonomous 
underwater vehicle designed to “provide researchers with a simple, low cost, rapid 
response capability” [16] to facilitate real-time oceanographic data collection in the 
marine environment. Although targeted at researchers, providing quick response to 
episodic events [16], the simple operational interface and relative low cost has significant 
appeal for military applications as well. Initial applications focused the development of 
REMUS vehicles to complement the Long Tenn Ecosystem Observatory (LEO-15), a 
cabled observatory system offshore from the Rutgers Marine Field Station near 
Tuckerton, New Jersey [17]. Multiple REMUS vehicles provide students and researchers 
continual access to underwater biological experimentation monitoring, ocean frontal 
system tracking, water column observation, sediment transport and bottom boundary 
layer growth and decay [16]. The capabilities have since expanded, paralleling the 
development of technology, and the REMUS is now used in a variety of missions, both 
civil and military, including [18]: 

• Hydrographic surveys 

• Mine Countermeasure (MCM) operations 

• Environmental Monitoring 

• Debris Field Mapping 

• Search and Salvage operations 

• Fishery operations 

• Scientific sampling and mapping 
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B. 


THE REMUS SYSTEM PACKAGE 


Available from Hydroid, Inc., the REMUS is provided in an all-inclusive package 
consisting of the REMUS vehicle, a ruggedized laptop with the REMUS GUI installed, 
power/data interface module, towfish transducer, a set of four navigation transponders 
and the REMUS Ranger. Figure 3 illustrates the REMUS 100 system package: 


B B B B 



Figure 3. REMUS 100 System Package, from [18] 


The design is a body-of-revolution, e.g. torpedo-shaped, based on a Myring hull 
profile [16], providing symmetry in the x-y and x-z planes. Table 1 provides the 
particulars for the original REMUS AUV (c. 2001): 
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Physical/Functional Area 

Characteristic 

Vehicle Diameter 

7.5 [in] 

Vehicle Length 

62 [in] 

Weight in air 

80 [lbs] 

External Ballast Weight 

2.2 [lbs] 

Operating Depth Range 

10 [ft] to 60 [ft] 

Transit Depth Limits 

328 [ft] 

Typical Search Area 

875 [yds] x 1093 [yds] 

Typical Transponder Range 

1640 [yds] 

Operational Temperature Range 

+32 [°F] to +100 [°F] 

Speed Range 

0.5 [knots] to 5.6 [knots] 

Maximum Operating Water Current 

2 [knots] 

Maximum Operating Sea State 

Sea Sate 2 

Battery 

1 [kW-hr] internally rechargeable Lithium- 
ion 

Endurance 

20 [hrs] at 3 [knots]; 9 [hrs] at 5 [knots] 


Table 1. Initial REMUS Functional and Physical Characteristics, from [13] 


The design of the vehicle is modular, greatly easing upgrade implementation 
based on new commercial technology as available. Due to this, the current iteration of 
the REMUS AUV in possession at NPS is significantly larger and more capable, 
presented in Table 2. 


11 




Physical/Functional Area 

Characteristic 

Vehicle Diameter 

7.5 [in] 

Vehicle Length 

107.25 [in] 

Weight in air 

145.60 [lbs] 

External Ballast Weight 

2.2 [lbs] 

Operating Depth Range 

10 [ft] to 60 [ft] 

Transit Depth Limits 

328 [ft] 

Typical Search Area 

875 [yds] x 1093 [yds] 

Typical Transponder Range 

1640 [yds] 

Operational Temperature Range 

+32 [°F] to +100 [°F] 

Speed Range 

0 [knots] to 5.6 [knots] 

Maximum Operating Water Current 

2 [knots] 

Maximum Operating Sea State 

Sea Sate 2 

Battery 

1 [kW-hr] internally rechargeable Lithium- 
ion 

Endurance 

20 [hrs] at 3 [knots]; 9 [hrs] at 5 [knots] 


Table 2. Current REMUS Functional and Physical Characteristics 


Of note, the length and weight of the vehicle vary significantly from its original 
counterpart. This is primarily due to the addition of the BlueView Forward Looking 
Sonar (FLS) package as well as the addition of the cross-body thruster package, both of 
which are shown in Figure 4. 
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Figure 4. BlueView Forward Looking Sonar Package 



Figure 5. REMUS AUV with Cross-Body Thruster Package and an uncovered 

BlueView FLS 

Other hardware and sensor upgrades available on the current NPS REMUS AUV 


are: 

• Inertial Navigation System 

• GPS Transponder 

• Acoustic Modem 
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Due to the near doubling of critical vehicle parameters, particularly length and 
mass, significant deviations of hydrodynamic and hull-form coefficients from those 
presented in [12] and [13] are expected. 

C. VEHICLE PROFILE 

Derivation of all necessary hydrodynamic coefficients requires a complete 
analysis of the REMUS vehicle, including hull profile, centers of mass and buoyancy, 
mass distribution, and control fin parameters. 

1. REMUS Hull Profile 

The REMUS vehicle is based on the Myring hull profile equations that describe a 
body contour with minimal drag coefficient for a given length-to-diameter ratio [12]. 
The REMUS vehicle profile is based on the following modified parameters defined with 
an origin at the vehicle nose: 

• a: the length of the nose section 

• b: the length of the constant-radius center section 

• c: the length of the tail section 

• d: the maximum diameter 

• n: an exponential parameter whose variation generates differing body 
shapes 

• 29: the included tail angle 

Nose shape is determined by the following modified semi-elliptical radius 
distribution [12]: 


i 



where r(E) is the radius as a function of axial position measured normal to the vehicle 
centerline, E is the axial position originating at the vehicle nose, a is the total nose length, 
and a offset is the missing length of the nose. 
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The tail shape is determined using the following equation [12]: 


•(H) = 


f i A 


W 


d- 


3 d tan# 


2c 2 


( e -/) 2 + 


d_ 

3 


where the forward body length, //is defined as: 


l f =a+b-a offset 



(6) 

(7) 


and, once more, r(S) is the hull radius as a function of axial position and S is the axial 
position. Figure 6 provides a relation of the parameters: 



Figure 6. Myring Profile relationships, from [12] 

Of note, the Myring profile includes a full tail length, while the REMUS vehicle 
itself has a snubbed tail, this difference is included as c 0 ff set . This hull profile does not 
include any external appendages such as the control fins, sonar transducer, external 
ballast and GPS transceiver. 

2. Centers of Mass and Buoyancy 

Due to the modularity of the vehicle, the mass of the REMUS can vary 
significantly between missions. For the purposes of this research, however, the mass of 
the REMUS vehicle is presumed constant. 


Parameter 

Value 

Units 

XCG 

-1.35 

[m] 

ycc 

0.00 

[m] 

ZCG 

-0.102 

[m] 


Table 3. Center of Gravity (datum at vehicle nose) 

15 





























Similarly, the vehicle buoyancy is constant with a slight positive buoyancy of 
approximately 0.5 [kg] as a failsafe to force the vehicle to surface in the event of a 
casualty. 


Parameter 

Value 

Units 

XCB 

-1.36 

[m] 

ycB 

0.00 

[m] 

ZCB 

0.00 

[m] 


Table 4. Center of Buoyancy (datum at vehicle nose) 


3. Mass Distribution and Inertia Tensor 


A critical assumption for the following hydrodynamic coefficient calculations 
involves determining the mass distribution and, in turn, the complete inertia tensor. As 
stated previously, the vehicle mass is unchanging. Traditionally, the mass distribution 
would sum the masses and locations of the individual vehicle components to provide a 
complete inertia tensor. As the vehicle components remain proprietary to their parent 
company, this work assumes a unifonn mass distribution based on a uniform mass 
density of the vehicle. The mass density is determined as: 

m 


Premus 


( 8 ) 


where m is the mass of the vehicle and V is the volume of the vehicle defined by: 

V = j hl nr (S) 2 dE (9) 

which is simply an integration of the instantaneous cross-sectional area of the vehicle 
over its length. Approximating the integral using Simpson’s Rule yields the following: 

j + 2^r(S 3 ) 2 +... 


(i/ 


V = 


/99 

3 


V 99 

I 

' n =0 


+2 


( nr ( s „-2 ) 2 ) + 4 () 2 ) + nr (S„ ) 2 


( 10 ) 
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where is the overall length sub-divided into 100 stations and ly—,,) is the vehicle 
radius at a given station. 


The inertia tensor is defined as: 



which, based on vehicle symmetry, can be reduced to the following: 

k, o o i 


(11) 


( 12 ) 


The product of inertia, I yz , is not neglected due to significant variations in vehicle shape 
in the forward and rear lengths. The general equations for remaining moments of inertia 
are: 


4 = j (4 + a) p A dv=Yjn n (>■; + z ]) (i3) 

y n =1 

4 = JV +z ' ) PA dV = Tj m n (4 + Z «) (14) 

y n =1 

4 = j(4 +/) p A dV = Yjn n (x„ 2 + y ]) (15) 

y n =1 

4 = j M pA dv =j( z >j p.A dv =4 = 2k k z «) (i6) 

y n =1 

Applying previous relations regarding vehicle radius and once more invoking axial 
symmetry, the general equations become: 


a. 


Moment of Inertia About the x-Axis 

4 = j (V + ) P.4 dV = 2k (4 + z « ) 


99 


M=1 


99 99 

=2k (4 s )! + 4 S )!)=2k ( 2r ( s )!) 

«=0 «=0 


(17) 
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where n is the station number, T{—) n is the hull radius at the specified station and m n is 
the mass of the station defined by: 


m n ~ Premus 



(18) 


b. Moment of Inertia About the y- and z- Axes 

Due to symmetry, the vehicle possesses identical moments of inertia about 
the y- and z- axes: 

I yy = \( xl+z 2 )PA dV = Z m " (*» + z «) 

V n =1 

4 = { (4 + y 2 ) PA dV = Z m n (4 + y 2 n ) ( 1 9) 

y n =1 

99 

4 = 4=Z«„( s « + r ( s )») 

«=0 

c. Product of Inertia Between the y- and z- Axes 

4 = JM pA dv =\{ z y) pA dv =4 = Z m » ) 

F 99 F " (20) 

4=4 = Z , M r ( s )«) 

n=0 


d. Resulting Moments of Inertia 


Parameter 

Value 

Units 

I XX 

1.154 


kg-m 


Iyy 

34.70 


kg ■ nf ” 


Izz 

34.70 


kg-m 2 ^ 


lyz 

0.5678 


kg ■ m 2 ” 



Table 5. Moments of Inertia (origin at vehicle half-length) 
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III. ELEMENTS OF THE GOVERNING EQUATIONS 


A. INTRODUCTION 

This chapter describes the assumptions and methods used to define the equations 
of motion, which model the motions experienced by the REMUS vehicle. 

B. MODELING ASSUMPTIONS 

The following assumptions apply to all aspects of this research: 

1. Vehicle Dynamics Assumptions 

• The REMUS AUVcan be modeled as a single rigid-body. In other words, 
the mass distribution of the vehicle remains constant and rigid-body 
dynamics and kinematics are valid. 

• The REMUS A UVpossesses constant mass. The vehicle neither gains nor 
loses any mass during the course of its operations. 

• Control response time is instantaneous. Change in control fin angle or 
thruster speed is immediate. 

• The REMUS AUV operates in a vehicle-centered coordinate frame. The 
coordinate frame originates at the vehicle’s center of buoyancy. 

2. Environmental Assumptions 

• The REMUS AUV is submerged in a “safe harbor” environment. That is, 
this model neglects wave effects and considers current effects as a linear 
additive. 

• Rotational effects of the earth are negligible. Acceleration components of 
gravity on the center of mass are constant. 

• Forces acting on the REMUS vehicle are gravitational and inertial. These 
forces are derived from hydrostatic, hydrodynamic and propulsive origins. 
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• Fluid properties are constant. That is, density (p), viscosity (//), kinematic 
viscosity (v), temperature and salinity are all constant for seawater at 
15.6°C. 

• Local variances in the earth’s gravitational field are ignored. 
Acceleration due to gravity is constant and assumed to be g = 9.81 



C. REFERENCE FRAME 

The following standard SNAME notation [19] presented in Table 6 will be used 
to describe the degrees-of-freedom, forces, velocities, positions and angles: 


Degree-of- 

Freedom 

(DOF) 

Motion Description 

Forces 
(F) and 
Moments 
(M) 

Finear 

and 

Angular 

Velocities 

Positions 

and 

Euler 

Angles 

1 

Surge - motions in the x-direction. 

X 

u 

X 

2 

Sway - motions in the y-direction. 

Y 

V 

y 

3 

Heave - motions in the z- 
direction. 

Z 

w 

z 

4 

Roll - rotations about the x-axis. 

K 

p 


5 

Pitch - rotations about the y-axis. 

M 

q 

e 

6 

Yaw - rotations about the z-axis. 

N 

r 



Table 6. Notation used for the REMUS AUV 


The first three DOFs describe the linear motion of the vehicle and the latter three 
describe the angular motion. It should be noted that, by convention, the positive x- 
direction is taken as forward, the positive v-direction is to the right, and the positive z- 
direction is down. This is commonly referred to as the North-East-Down [NED] 
coordinate frame. Furthennore, it is convenient to annotate the quantities in Table 6 into 
the following vectors: 
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(21) 


t = [X,Y,Z,K,M,N] T 
v = [u,v,w,p,q,r] T 
x = \x,y,z,(/>,9,y/] T 

To properly describe the motions of the REMUS vehicle, two orthogonal 
coordinate frames are required. The first is a global coordinate frame, usually considered 
to be earth-fixed. The second is a body-fixed coordinate frame, attached to the REMUS 
vehicle itself at the center of buoyancy (CB) and oriented to NED. Figure 7 shows the 
relation between the two coordinate frames and their respective forces, moments, 
velocities, and positions: 



BODY-FIXED 
COORDINATE SYSTEM 


INERTIAL OR EARTH-FIXED 
COORDINATE SYSTEM 


SWAY (2): v,' 
PITCH (5): q, 



HEAVE (3): w,Z 


SURGE (1):u,X 


. 

YAW (6): r.N 


ROLL (4): p,K. 



Uuucu NT 


ghfgjfe 'Vr'; *1 • 



Figure 7. Coordinate Frame Relations, from [12] 


Utilization of Euler Angles enables the transformation from one coordinate frame 
to another in the Cartesian coordinate system. The use of Euler angles requires each 
coordinate axis to be transformed individually. That is, the body must first be 
transfonned about the z-axis, then transformed about the y-axis, and finally about the x- 
axis. The individual rotation matrices are: 
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cos(^) 

sin(^) 0 

K v = 

-sin(^) 

cos(^) 0 


0 


0 Ij 


cos(#) 

0 

-sin (#) 

Ke = 

0 

1 

0 


sin ( 9 ) 

0 

cos(#) 


[1 0 


0 

Ke = 

0 cos(^) 

sin(^S) 


0 sin(^) 

cos(^S) 


The total rotation matrix is the multiplication of the three: 


(22) 


(23) 


(24) 


R{(f>,9,y/) = 

cos (y) cos (#) -sin(^)cos(^) + cos(^)sin(#)sin(^) sh^^cos^ + cos^^os^sin^) 
sin (i//) cos (9) cos(^)cos(y) + sin(^)sin(#)sin(y) -cos(^ r )sin(^) + sin(^ r )sin(6 , )cos(^) 

-sin(#) cos (9) sin (0) cos ( 9) cos (</>) 


u 


~X 

V 

= R(</>,9,if/) 

Y 

w 


Z 


The above rotation matrix only enables the transformation of the body velocity vector to 
the global velocity vector; transformation of the angular velocities must be done as 
follows: 


P 




~ 0 ~ 


" 0 " 

9 

- R ,„ 

0 

+ R xrf R y,0 

9 

+ R x^ R y,0 R z,y, 

0 

r 


0 


0 


y_ 


p 


"l 

0 -sin(0) 


9 

= 

0 

cos(^) sin {(/)) cos [9) 

9 

r 


0 

-sin(^) cos [(/)) cos (9) 

¥ 


(26) 


Multiplying the angular velocity vector by the inverse combined rotation matrix 
yields angular accelerations: 
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(j) 1 sin(y)tan($) cos(y)tan(<9) p 

6 = 0 cos(^) -sin(^) q (27) 

ys sin(^) cos(^) r 

cos(6*) cos(6*) 

It is important to recognize the singularity at pitch angles of (?=±90°. Generally 
speaking, underwater vehicles do not usually operate near these singularity points. The 
REMUS AUV, in particular, is designed to operate outside these singularities. 
Furthermore, the application of the cross-body thrusters provides any required direction 
change through both linear or angular forces and moments, instead of relying solely on 
angular moments applied by fin control surfaces. 


EQUATIONS OF MOTION 

The equations of motion are an extension of Newtonian physics: 

yp = m-d 
Y j M=I-g) 


where 



(29) 


(30) 


where M RB and M A are the rigid-body and added mass matrices, respectively and M is the 
complete mass matrix. 


1. Vehicle Dynamics 

Recalling the assumptions made in Chapter III.B, the above vector representation 
is broken down into its respective degrees-of-freedom below: 
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X= m\u -vr+wq-x CG (q 2 +r 2 ) +_y CG (pq-r) +z CG ( pr + q) 
Y=m\v-wp+w-y CG (r 2 +p 2 )+z CG (qr-p)+x CG (qp+r) 
Z = m \w-iiq+vp-z CG (p 2 +q 2 )+x CG (ip-q)+y CG (rq+p) 

K = I xxP + (L - I yy)qr-{r + pq)I„ +(r 2 -q 2 )l yi + (pr-q)l x 
+ m \_y CG ( w-uq + vp ) - z GG (v - wp + ur )] 

M = I„q + (4 - 4 )rp-(P + c l r ) 4 +(p 2 -r 2 )l a +{qp-r)l x 

+m[z CG (w -vr + wq)-x CG (w-uq + vp)^ 

N = IJ + (l yy -I xx )pq~{q + rp)l yz + (q 2 -p 2 ) 4 +(rq-p)l 2 
+m [x CG (v - wp + wr) - _y CG (li - vr + wp)] 
where in is the vehicle mass and /,, is the following complete inertia tensor: 


1 = 


I 

xx 

-I. 


-I -I 

xy x: 


-I -I 


(31) 

(32) 

(33) 

(34) 

(35) 

(36) 


(37) 


The first three equations model the translational motions of the vehicle and the second 
three model the rotational motions. Recall that the origin of the body-fixed coordinate 


frame is located at the center of buoyancy, CB= [x c/; y CB z CB J ; as such, the equations 

neglect these zero-valued terms. Furthermore, recalling the assumption that the REMUS 
vehicle possesses two planes (xy- and yz-) of symmetry, the products of inertia I xy , I yx , I xz , 
and I zx are negligible, reducing the inertia tensor to: 


/ = 


0 

0 





(38) 


This reduces the equations of motion to: 


X=m 


Y —m 


Z — m 


u-vr+wq -x CG (q 2 +r 2 ) +y CG (pq-r) +z CG (pr + q) 
v-wp+w-y CG (r 2 +p 1 )+z CG (qr-p)+x CG (qp+r ) 
w-nq+yp-z CG (p 2 +q 2 )+x CG (ip-q)+y CG (rq+p) 


(39) 

(40) 

(41) 
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K = I xx p + (4 - 4 ) qr + (r 2 - q 1 ) I yz + m[y CG ( w-uq + vp)~z CG (v -wp + urj\ (42) 
M = I y y q + { I xx - I zz) r P + {<lP-r)lv +m[z CG {u-vr + wq)-x CG (w-uq + vp)] (43) 
N = IJ + (4 -4 )pq~(q + rp)l yz + m\x CG (v -wp + ur ) - y CG (u -vr + wq)\ (44) 


Based on vehicle symmetry, ballasting characteristics, and internal layout, variances in 
the location of ycG arc negligible. This further reduces the equations of motion to: 

X=m\ 


Z—m 


w- 


+r 2 ) 

+z CG (pr+q) 



(45) 

~P) 

+ x CG (qp + r)] 



(46) 

1 2) 

+ * , 

+ 

a* 

-3 

1 

1 _ 1 



(47) 

4- 

m[z CG (v-wp 

+ ur)] 


(48) 

CG {P 

-vr + wq)-x CG ( w-uq + vp )] 

(49) 

[ y Z +m[x CG (v-wp+ur)] 


(50) 


2. Vehicle Mechanics 

In the preceding equations of motion, the external forces and moments, described 
in Equations (45) through (50) are described in terms of vehicle coefficients. These 
coefficients are a combination of theoretical equations and empirical formulations. The 
actual coefficient values are derived in the following chapter. 
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IV. COEFFICIENT DERIVATION 


A. INTRODUCTION 

This chapter derives the coefficients representing external forces and moments 
upon the REMUS vehicle. While the standardized reference frame [19] remains in use, 
the correlating coefficients will retain their respective dimensions to facilitate comparison 
between this work and Prestero [12], contrary to the standard conventions presented in 
[19]. 

B. ADDED MASS 

A vehicle’s added mass represents the mass of the water column opposing the 
vehicle’s line of motion, and is greatly influenced by the shape of the vehicle. Typically, 
the REMUS vehicle normally travels along a longitudinal axis of motion, introducing 
added mass to the nose and body from surface friction. Its slender profile, however, 
yields a relatively small added-mass tensor in this direction. Application of the cross¬ 
body thrusters, however, induces a motion subjecting the vehicle’s entire length to a more 
resistive added mass tensor. Fossen [20] provides the full added mass matrix for an 
asymmetric submerged body: 


X, 

X, 

x* 

X P 

X i 

X f - 

Y* 

Y, 

Y * 

Y P 

Y P 

Y > 

z„ 


Z* 

Z P 


Z< 

K u 


K * 

K P 

K * 

K> 

K; 

M v 


M P 


M, 

V, 

N v 

N* 

N P 


Nr 


Because the REMUS vehicle is symmetric in the x-y (top-bottom) and x-z (port- 
starboard) planes, the added mass matrix reduces to: 
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** 

0 

0 

0 

0 

0 

0 

Y* 

0 

0 

0 

Y, 

0 

0 

Z+ 

0 


0 

0 

0 

0 

K. 

P 

0 

0 

0 

0 

M. 

vv 

0 

M. 

q 

0 

0 


0 

0 

0 



(52) 


Expansion of the vector equation for fluid kinetic energy, T A — 


\?M A v, yields: 


2 T A = X ll u 2 + v(Y r v+Y r r) + w(Z„w+ 


Zjl) + K pP 2 + 9 (Mjv+Mfl) +r(N ( v + N f r) (53) 


Considering Kirchoffs Equations in component form [20]: 


d 8T. 8T, 


dt du 
d dT A 

dt dv 
d ST, 


= r- 


q- 


dT s 


dv dw 
dT. dT, 


- — r- 


dt dw 


= q- 


dw du 
dT. dT, 


du 


d dT A 
dt dp 


■ = w- 


8T, dT. 


d dT 


dt dq 
d dT A 
dt dr 


A _ 


dv 

dT\ 

dw 


-v- 


dw 

dT\ 

du 

dT, dT, 


- + r- 


-w- 


+ p 


= v- 


du 


— u - 


dv 


+q 


dv 

El. 

dq 

El 

dr 

dT, 


-X A 


-Y, 


Z, 


-r- 


,E 

dr 

El 

dp 
dT. 


-K a 


-M 


dp dq 


-N a 


(54) 

(55) 

(56) 

(57) 

(58) 

(59) 


Substitution of Equations (54) through (59) into Equation (53) yields the complete set of 
equations for added mass: 


Ya 

Za 

K a 

M a 

N a 


X u u 2 + Z^wq + Z q q 2 - Y v vr - Y r r 2 
Y v v + Y r r + X.ur-Z.}vp -Z pP q 
+ Z q q - x i,m + Y,vp + Y.rp 

K pp 

M n W + M.q - (Z w - X, ) uw- Y.vp + (k.~ N,)rp - Z.uq 
NX + N,r -{X.- Y.)uv + Z.wp- (X p - M .) pq + Y r ur 


(60) 


The remainder of this section provides a breakdown of the added mass coefficients. 
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1 . 


Axial Added Mass 


The calculation for axial added mass approximates the REMUS hull-form as an 
ellipse for which half the vehicle length, /, is the major axis and half the vehicle’s 
maximum diameter, d max , is the minor axis. Blevins [21] presents the following 
empirical formula for the axial added mass of an ellipsoid: 


x,= — 


4 apn 

rn 

ro 

3 


l 2 ) 


(61) 


or 


X =■ 


Aftpn 




v ^ J 


(62) 


where p is the local fluid density and a and [> are empirical parameters measured as a 


function of the vehicle length to diameter, l/d max , ratio. 


ud 

a 

P 

0.01 

- 

0.6348 

0.1 

6.148 

0.6148 

0.2 

3.008 

0.6016 

0.4 

1.428 

0.5712 

0.6 

0.9078 

0.5447 

0.8 

0.6514 

0.5211 

1.0 

0.5000 

0.5000 

1.5 

0.3038 

0.4557 

2.0 

0.2100 

0.4200 

2.5 

0.1563 

0.3908 

3.0 

0.1220 

0.3660 

5.0 

0.05912 

0.2956 

7.0 

0.03585 

0.2510 

10.0 

0.02071 

0.2071 

00 

0 

- 


Table 7. Added Mass Parameters a and (3, from [21] 
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Figure 8 plots a and /? as a function of lld ma x- 


Blevins Parameters 





6.0000 

■ 



♦ Alpha 

■ Beta 


- Power Regression - Alpha 

- Exponential Regression - Beta 


1 


\ 

1 

-— i — . 

0.0000 • w ♦ i > 

0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00 

I/d 


Figure 8. Blevins Parameters 


/ 

Since ~ = 14.34 a which lies well outside the boundaries provided by Blevins, 

inax 

regression lines were plotted to extrapolate the required values for a and p yielding the 
following relations: 

,1.235 


f -U = 0.4466 

( l \ 

\ «max ) 

V^max y 


^(^r) = °- 5720831e 

Appendix B tabulates all final coefficient values. 


-0.1148249 


(63) 


2. Cross-Flow Added Mass 

Cross-flow added mass is computed using strip theory, a numerical integration 
technique. This term requires the combination of both cylindrical and cruciform hull 
cross sections. Blevins [21] indicates that the added mass per unit length of a cylindrical 
slice is: 
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( 64 ) 


m a (E) = 7rpr(E) 2 

where p is the density of the surrounding fluid and r(E) defines the hull radius as a 


function of axial position. Blevins also defines the added mass of a circle with symmetric 
fins as: 


m 


of 


(S) = npa 


\ '•(s)V(Er 


V 


a 


a 


(65) 


where ay is the maximum fin height measured from the centerline of the vehicle. 


Integrating the above two equations over the length of the vehicle yields the 
following equations for added mass: 


K = -£' m a (E)dE-f* m af (E)dE- Q m a (S )dE 

(66) 


(67) 

z w =l 

(68) 

II 

nT 

(69) 

M w = £' E-m a {E)dE-^E-m a (E)dE-^ E-m a (E)dE 

X t * x f * x f2 

(70) 

' / Vl 

(71) 


(72) 

A 

II 

(73) 


Appendix B describes and tabulates all final coefficient values. Table 13 in Appendix A 
provides the limits of integration. 


3. Rolling Added Mass 

The derivation of the rolling added mass term assumes that the cylindrical shape 
of the hull and any small protrusions thereof have no effect on the rolling added mass 
tenn. Based on this assumption, the only section of the vehicle contributing any 
appreciable rolling added mass is the tail section containing the vehicle’s control surfaces 
(fins). Blevins [21] offers the following empirical formula for added mass of a rolling 
circle with multiple equally spaced fins: 


31 



( 74 ) 


V f X/2 2 „ 4 7P—' 

a . = —pa di .i 

J v n 

where a is the average fin height above the centerline of the vehicle. Table 13 and Table 
14 in Appendix A provide the limits of integration and fin parameters, respectively. 

4. Added Mass Cross-terms 

The remaining added mass terms resulting from added mass coupling can be 
determined from calculations previously derived following the relationships presented in 
Table 8. The added mass cross-tenns are grouped by row according to the respective 
degree of freedom and represent added mass along a degree of freedom due to coupling 
from off-axis velocities, both linear and angular. 


N 

II 

^7 

T =Z. 

II 

II 

k 

wq w 

qq q 

W V 

rr r 

* 

II 

7 

II 

i 

ii 


ur u 

wp W 

pq q 


z =-x. 

>- 

II 

N 

Z =Y. 


uq u 

vp V 

rp r 
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Table 8. Added Mass Cross-tenn Relations 


The values for the added mass cross-terms are presented in Appendix B. 

C. DRAG 

Drag is a resistive force upon the vehicle that opposes the vehicle’s line of 
motion, and is heavily influenced on the shape of the body, particularly the exposed area 
along the axis of motion. The governing equation for fluid drag is [22]: 

F D =\pC D AV 2 (75) 

where F D is the viscous drag force, p is the density of the surrounding fluid, C'd is the 

drag coefficient of the body, A is the area exposed along the axis motion, and V is the 
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body’s velocity along the axis of motion. Since the viscous drag force always opposes 
the line of motion, it is necessary to consider the squared velocity tenn, V 2 , as V|V|, 
resulting the following governing relation: 

Fd=\pC d AV\V\ (76) 

yielding the appropriate sign. 


1. Axial Drag 

Axial drag is expressed by the following empirical relationship: 


* D =~ 


1 




pC D A F 
v2 j 


u\u\ 


yielding the following non-linear axial drag coefficient: 

A =--pC n A F 

U\U\ D t 


(77) 


(78) 


meaning 


X D =X^\u\ (79) 

where p is the density of the surrounding fluid, Cd is the vehicle’s axial drag coefficient 
and Af is the vehicle’s frontal area. 


Hoerner [23] presents the following empirical fonnula to calculate the axial drag 
coefficient, Cd- 



(80) 


where the skin-friction drag coefficient, C/, is determined using an approximation to the 
Schoenherr Equation [23]: 


or, more clearly: 



3.461og(Re)-5.6 


C, 


^3.46 log (Re)-5.6 J 


( 81 ) 


(82) 
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Hoemer states that the Schoenherr Equation is valid throughout the turbulent 
region of Reynolds numbers [23]: 

_ VI 

Re = — (83) 

v 

where V is the vehicle’s velocity, / is the length, and v is the kinematic viscosity of the 
surrounding fluid. Figure 9 shows the range of Reynolds numbers for the REMUS 


vehicle in seawater at 15.6°C at speeds of: 0.25 


in/ 


< V< 2.57 




The range of Reynolds numbers particular to the REMUS vehicle clearly fall 
within the standard values of the turbulent flow region, indicating the validity of using 
the approximation to the Schoenherr equation. Figure 10 and Figure 11 illustrate the 
trends for the drag coefficient, Cd, and the axial drag coefficient, X u \ u \, relative to the 
standard velocity range shown in Table 2. 
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Figure 10. Drag Coefficient, Cd 



Figure 11. Axial Drag Coefficient, X,, 
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The data trends follow an intuitive physical meaning - as velocity of the vehicle 

approaches 0 , Cd approaches +oo and X u \ u \ approaches -oo, consistent with the notion 

that a body at rest experiences an infinite drag resistance opposing its motion. 
Conversely, as the speed of the vehicle approaches infinity, the drag coefficients 
approach an asymptotic value. 

Final values for both the skin-friction drag coefficient, C/, the drag coefficient, 
C D , and the axial drag coefficient, X u \ u \, are provided in Appendix B. 

2. Cross-Flow Drag 

Cross-flow drag is the summation of drag forces normal to the fluid flow; in the 
case of the REMUS vehicle, this includes both hull- and fin-induced cross-flow drag. As 
with axial drag, the shape of the body heavily influences the amount of cross-flow drag 
experienced. The method used to calculate cross-flow drag is analogous to the method 
used in Equations (66) through (73), however the accuracy when applied to viscous terms 
can be off by as much as 100% [12]. Inclusion of these terms, however, yields a fully 
defined model whose coefficients can be corrected with experimental data and 
observations at sea. 

The non-linear cross-flow drag coefficients are calculated as follows: 

K ,H = Z H = ■ - \pc* 2r 0 3 )' d ' ! =- 2 ( \pS f ,c« ] (84) 

m -H = ~ N M = \ 2s ''( s ) ds - 2 *», () pSf„Cjf ) (85) 

2 V2 J 

Y M-~ z q \ q \ = -^P c &{ i rA22s | s k( s )^ s - 2 ^„|^„|f^P^A/] (86) 

M 11 Z yx t \Z J 

"*|= (87) 

where p is the density of the surrounding fluid, c dc is the drag coefficient of a cylinder, 

S is the axial position (shown in Figure 6 r(S)is the radius as a function of axial 

position (calculated in Equations (5) and (6)), x fm is the axial fin-post location, S fm is the 

fin planar area and c df is the cross-flow drag coefficient of the fins. 
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The cylindrical cross-flow drag coefficient, c dc , is estimated by [Hoemer] to be 
1.1. The fin cross-flow drag coefficient, c df , is derived using the following formula 

developed by Whicker and Fehlner [24]: 

c df = 0.1 + 0.7t (88) 

where t is the fin taper ratio. 

Appendix B presents the final non-linear force and moment coefficients and Table 
14 presents the REMUS fin parameters. 


3. Rolling Drag 


The rolling resistance of the REMUS vehicle will be approximated using the 
assumption that the primary contribution comes from the control fins, ignoring any other 
hull protrusions. This implies that: 

F n =(y a p\p\ (89) 

D roii \ VV/ avg J avgy \r\ V / 

where F D is the drag force due to roll, Y is the fin component of the vehicle cross-flow 

drag coefficient calculated by reducing the limits of integration from Equation (84) 
include only the hull section possessing fins, and r mean is the mean fin height above the 


vehicle centerline. This defines the rolling drag coefficient as: 


K =7 a 3 

p\p\ w y avg 


(90) 


This method provides only a best estimate of rolling drag; experimental data would 
provide a better result. 


D. LIFT 

As a vehicle moves through a fluid at a specific angle of attack, it will experience 
lift as the flow separates around the vehicle and a pressure drop occurs. Depending on 
the angle of attack of the vehicle, this pressure drop will induce a pitching moment about 
the center of buoyancy and either force the vehicle to submerge deeper or begin to 
surface. In underwater vehicles, body lift, the lift force generated by the angle of attack 
and shape of the hull itself, is significant compared to the forces generated by the control 
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surfaces. Specifically, the body has a much larger area exposed to the flow. In aerial 
vehicles, on the other hand, the control surfaces are much larger, reducing the impact 
body lift has on the dynamics of the vehicle. 


1. Body Lift Force 


The body lift force is governed by the standard lift equation [25]: 

F L =^pC L AV 2 


(91) 


where F L is the viscous lift force, p is the density of the surrounding fluid, Cl is the lift 
coefficient of the body, A is the area exposed along the axis of motion, and V is the 
body’s velocity along the axis of motion. Nahon [26] proposes the following linear 
relationship between body lift coefficient and the vehicle’s angle of attack: 


„ _2(k 2 -k x )S» 


(92) 


where (k 2 -ki) is an added mass factor detennined by Munk [27], So is the cross-sectional 
area of the body station where the flow ceases to be potential and V is the hull volume. 
Munk’s added mass factor, (k 2 -ki), is a difference between transverse, k 2 , and 
longitudinal, ki, added masses for elongated surfaces of revolution. Table 9 provides 
Munk’s data: 


l/d 

ki (Longitudinal) 

k 2 (Transverse) 

k 2 - ki 

1.00 

0.500 

0.500 

0.000 

1.50 

0.305 

0.621 

0.316 

2.00 

0.209 

0.702 

0.493 

2.51 

0.156 

0.763 

0.607 

2.99 

0.122 

0.803 

0.681 

3.99 

0.082 

0.860 

0.778 

4.99 

0.059 

0.895 

0.836 

6.01 

0.045 

0.918 

0.873 

6.97 

0.036 

0.933 

0.897 

8.01 

0.029 

0.945 

0.916 

9.02 

0.024 

0.954 

0.930 

9.97 

0.021 

0.960 

0.939 

00 

0.000 

1.000 

1.000 


Table 9. Munk’s Added Mass Coefficients, from [27] 
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Use of Nahon’s relation requires two assumptions based on the dimensions of the 
REMUS vehicle: 

(*2-*0 = 1-000 (93) 

$)=<L (94) 

The first assumption is based on the fact that the REMUS vehicle has a considerably 
larger l/d max ratio falling well outside Munk’s empirical values. The second assumption 
dictates that the flow around the REMUS vehicle is steady through the initial contact with 
the vehicle’s nose. Considering the smoothness of the hull-form at the nose, the 
uniformity of the diameter through the majority of the vehicle’s length and the locations 
of various appendages, this assumption is reasonable. This implies the following 
relationship for the overall body lift force: 

< 95 > 

meaning 

1 (2d \ 

Y u V , = Z U,1 =-P A (96) 

z v v /3 j 

Appendix B tabulates all final coefficient values. 

2. Body Lift Moment 

Hoerner estimates the center of the lift force for an ellipsoid is “in the vicinity of 
30% of the chord” [28]; in this case the chord is the overall length of the vehicle. This 
defines the center of pressure as: 

x q, =-0-70/ (97) 

meaning that the body lift moment is: 

K, =#«* (98) 

Appendix B tabulates all final coefficient values. 
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3. 


Fin Lift Force 


Attitude control, specifically rotations about the axes (i.e. yaw, pitch, and roll), 
are controlled in part by a pair each of rudders and stern planes arranged in a cruciform 
pattern at the aft end of the vehicle. Each control pair is slaved, that is individual rudders 
are not individually controlled. For vehicle control fins, the empirical formula for fin lift 
is given as [12]: 



(99) 


where F L is the lift force from the fin, C L is the fin lift force coefficient, S f m is the fin 
planform area, 8 e is the effective fin angle in radians, and V e is the effective fin velocity. 


The fin lift coefficient, C L , is a function of the effective fin angle of attack, a. Hoerner 
[28] provides the following relation for the fin lift coefficient as a function of a: 




dC, 


da 


1 1 

2 a 7i 7i(AR e ^ 


( 100 ) 


where CL is a form factor for foil sections in undisturbed flow found by Hoemer to be “in 
the order of 0.9 [28]” and AR e is the effective aspect ratio of the fin, given by the formula 
[ 12 ]: 


AR=2(AR) = 2 


( 


V 


J 


( 101 ) 


where bf m is the fm span. 


Since the location of the control fins does not coincide directly with the origin of 
the vehicle coordinate system, the fins experience a slightly different effective velocity 
[ 12 ]: 


U fin = ll + Z fi n Cf~yfin r 

(102) 

Vfin= V+X fin r - Z finP 

(103) 

W fin= W+ y fi nP- X fin ( l 

(104) 


where Uf m , Vf in , and Wf m are the effective fm velocities and Xf m , yf m , and Zf m are center-post 

locations of the control fins in the body-referenced coordinate system. Since the control 
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fins are placed in a standard cruciform pattern and because translational terms dominate, 
terms involving yj}„ and Zf m are negligible. 

The very nature of the control surfaces dictates that they experience an effective 
angle of attack relative to their orientation to the vehicle and the vehicle’s orientation to 
the fluid. 



Figure 12. Effective Rudder Angle of Attack, from [12] 



Figure 13. Effective Stern Plane Angle of Attack, from [12] 


These effective angles of attack are expressed as [12]: 

S ,,,= 5 r-Pr, ( 105 ) 

K= S s+P,„ (106) 

where d r and d s are the control fin angles in the body-referenced coordinate system, 
fi rf{ and P S ] are the effective angles of attack of the fin zero plane as shown in Figure 12 

and Figure 13. Assuming small angles and carrying the assumptions made for relative 
fin velocities, these effective angles are expressed as [12]: 

( v+ v) 


/L = 


r eff 


V 

jin 


t j- 

jin 


(107) 
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( 108 ) 


g _ w fi« A w ~ x fl» q ) 

•" u fm ~ u 

Substitution of Equations (100) through (108) into Equation (99) yields the following 
relationships governing fin lift: 


K =^pC Lfi S fl „ \}f8 r -u\’-x flll ur\ 

(109) 

1 r •> 

Z = --pC Lfi S fm \u8 s +uw-x fm uq\ 

(110) 

meaning that the fin lift coefficients are: 


Y uuSr = ^P C L fi fifln 

(HI) 

Z , =-— pC, S f 

UUOS ^ y tW 

(112) 

Y =-Z =- — pC, S r X r 

ur f uwj■ ^ ' Lf in fin fin 

(113) 


Appendix B tabulates all final coefficient values. 

4. Fin Lift Moment 

The fin lift moment is a simple extension of the fin lift force over the distance of 
the control fins from the body-referenced origin: 


M, =F, x, 

L jm L Jm A" 

(114) 

This results in the following equations for fin lift moments: 

1 r ? i 

K =-pC Lfi S fln \u 8 s +uw-x fln uq\x fin 

(115) 

Y r =\p C I, m S fin [ l{2S r ~ UV ~ X fi n Ur ] X /, n 

(116) 

further resulting in the following fin lift moment coefficients: 

M-uuds 2 P^L )tn ^'fin X fin 

(117) 

N U uSr = ~ P^L fi „ Sfin X fin 

(118) 
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(119) 


^uq f N w. f 2 PC Lfl S fin X fin 

Appendix B tabulates all final coefficient values. 

E. HYDROSTATICS 

The effects of weight and buoyancy impart the observed hydrostatic forces on the 
vehicle. Assuming m is the mass of the vehicle; requiring weight, W = mg, where g is the 
local acceleration due to gravity. Buoyancy, the restoring counter to weight is expressed 
as: 

B = pgV (120) 

where p is the local fluid density and V is the volume of the vehicle. The REMUS AUV 
is designed to be slightly positively buoyant; a safety feature that guarantees eventual rise 
to the surface in the event of a malfunction. The combined hydrostatic forces and 
moments can be expressed as: 

^m = fcG~fcB (121) 

^hs ~{ r CG X fco)~ ( r cs x ^ce) (122) 

Applying the Euler angle transformations matrix, derived in Equation (25) 

from Chapter III yields: 

“ 0 " 

fco=K, 0 0 23 ) 

W 

“ 0 “ 

fcB=Kv ° ( 124 ) 

B 

Expansion of this equation set yields the following set of hydrostatic forces and moments 
on the vehicle: 
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( 125 ) 




'X m 


-(W-B) sintf 



Yhs 


(W -5) cos 0 sin (j) 

F hs ' 


^HS 


(W - B)cos6cos(j) 

F^HS _ 


K hs 


-(v cc W - y CB B) cos 0 cos (/> - ( z CG W - z cb B ) cos 6 sin (f) 



M hs 


- (z CG W - z cb B ) sin 9 - (x CG FF - x cb B) cos 0 cos (f) 



N hs_ 


~ (*CG W - x cb B ) cos 6 sin (f) - {y cc W - v cb B )sin 0 


F. PROPULSION MODEL 

Propulsion for the REMUS vehicle is provided from 5 thrusters: a 3-bladed main 
propeller of diameter D prop = 0.1143 [ml, and four cross-body thrusters of diameter 

D t hn,ster = 0.0381 [m]. Modeling thrust for the REMUS vehicle requires a relationship 

between the input parameter of RPM and the output parameter of measured thrust. No 
governing empirical or analytical relationship between the two exists since thrust is a 
function several variables [29]: 

F r = f {p, D,V a , g,n, p, fi) (126) 

where 

F t = Thruster force 
p = Density of the fluid 
D = Diameter of the propeller 
V A = Speed of advance of the propeller 
g = Acceleration due to gravity 
n = Speed of rotation 
p = Fluid pressure 
p = Viscosity of the fluid 

Recalling the assumptions made in Chapter III, the relationship becomes simply: 

F r =f(n) (127) 

This is consistent with the control interface for the cross-body thrusters in the included 
REMUS GUI. The empirical relationship was detennined analytically in a tow tank test 
at the Naval Postgraduate School. 
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1 . 


Thruster Force Measurement 



Suspending the REMUS vehicle from a crane and using a FUTEK USB Strain 
Gage, shown in Figure 14 provided accurate measurements of the forces imparted by the 
vertical thrusters. 


Figure 14. FUTEK USB Strain Gage 

The two thrusters were controlled in increments of 150 RPM or 3% of maximum 
RPM from 0 RPM to 5000 RPM. It is important to note that the pairs of thrusters are not 
slaved; that is, incremental change to one does not affect the other. This resulted in a 
slight lag as settings were altered. Measurements were collected after a lengthy period of 
stabilization to mitigate the effects of the un-slaved system and errors therein. Figure 15 
plots the measured thrust as a function of rotational speed. 
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Measured Thruster Data 



0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 

Thruster Speed (RPM] 


-Measured Thrxister Data 


Figure 15. Measured Thruster Data 


The measured data adheres closely to the following 4 th order polynomial regression: 

F t = -(8.0712-lO" 14 • n 4 ) +(5.6674-lO -10 • n 3 )-(4.6335-10“ 7 -n 2 ) +(7.483 M0“ 4 •«) 

(128) 

-(1.0795-10" 2 ) 

which is used in the REMUS model to provide an accurate assessment of thruster power 
from rotation speed input. This relation is assumed equivalent for all four cross-body 
thrusters, as all four thrusters are of identical size. Furthermore, since a polynomial 
regression is unbounded, the 4 th order regression is only valid between the minimum and 
maximum rotation speeds. 


2, Thruster Torque Calculation 

Thruster torque is calculated using the following relation: 

T=F 7 jc Ioc (129) 
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where T is the calculated torque, F T is the measured thrust, and xi oc is the location of each 
thruster relative to the vehicle’s center of buoyancy. Plotting torque as a function of 
rotational speed results in Figure 16. 



The minor variations in torque are due to the differing location of each thruster 
relative to the center of buoyancy. Similar to the thruster data, the torque calculations 
adhere to a 4 th order polynomial fit, resulting in the following set of equations: 


T flt =-(3.128-10“ 14 -n 4 ) + (2.196-lO' 10 V)-(l.795•10" 7 -n 2 ) + (2.900-10” 4 -n) 
-(4.183 -10 -3 ) 

T fvt =-(2.954 • 1 O' 14 V) +(2.074-lO' 10 V)-(l.696-l(T 7 -n) + (2.73910 4 -n) 
- (3.951 • 10 -3 ) 

T slt = ( 3 .196 • 10“ 14 • n ) -(2.244 ■ 10" 10 • n ) + (l.835 • 10' 7 • n 1 ) -(2.963 ■ 1(T 4 • n ) 

+ (4.275 • 10" 3 ) 


(130) 

(131) 

(132) 
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T SVT = (3.350■ 10 -14 • » 4 )-(2.352-lO" 10 •« 3 ) + (l.923'10“ 7 -/7 2 )-(3.106■ 10 -4 -n 
+ (4.480 ■ 10" 3 ) 

Again, these equations are only valid between the minimum and maximum rotation 
speeds. 



G. COMPLETE HYDRODYNAMIC TERMS 

Description and tabulation of all non-linear force and moment coefficients derived 
in this chapter is provided in Appendix B. 
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V. COMPLETE MODEL AND TESTING 


A. INTRODUCTION 

This chapter combines the equations of motion discussed in Chapter III with the 
hydrodynamic coefficients derived in Chapter IV to provide a complete, six degree-of- 
freedom model for the REMUS vehicle. The SIMULINK model is adapted from the 
MATLAB model of the NPS Phoenix AUV provided in the Marine Systems Simulator 
Guidance, Navigation and Control (MSS GNC) toolkit [30], 

B. COMBINED NON-LINEAR EQUATIONS OF MOTION 

The combination of equations for rigid-body vehicle dynamics and the equations 
involving external forces and moments provide the complete set of non-linear equations 
governing the motion of the REMUS vehicle. These non-linear equations incorporate all 


the hydrodynamic 

coefficients and 

map 

the control 

input 

vector, 

u = [S r S p n FLT n 

T 

SLT ^ FVT ^ SVT ^ ? 

a series 

; of applied forces on the vehicle. 

The equations 

are presented 

below 

following 

the 

fonn 

of V F = y F 

body external 

(or Ys M body =Ts M external 

) = 





Surge 


m 


■ vr + wq - x CG ( q 2 3 + r 2 ) + y CG (pq - r) + z CG (pr + q ) 


— Xf,c H" X | \U \u + X u 

no u\u\ i i u ^ 234 ^ 


+X wq + X q 2 +X vr + X r 1 + X 

wq j qq 1 vr rr prop 


2. Sway 

v-wp + ur - y cc ( r 2 + p 2 ) + z CG (qr -p) + x CG (qp + r ) 
+Y-v + Y.r + Y ur + Y wp + Y pq + Y uv + Y ,u 2 S +Y 

v r ur wp r pq* 1 uv uuo,. r prop 


= Y HS +7, ,v v +7 ,r r 

v|v| | | r\r\ | 


m 


3. Heave 

w-uq + v p- z CG (p 2 + q 2 ) + x CG (rp -q) + y CG ( rq + p) = Z HS + Z , yv\w\ + Z^q\q\ 


+Z w w + Z.q + Zuw + Zuq + Zvp + Zrp + Z uug u 2 S v + Z c 


(136) 
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4. 


Roll 


J xx p + (4 - 4) V + 4 ('" 2 - q 2 ) + m [ y CG ( w-uq + vp ) - z CG (v - wp + Mr)] = K m 


+ 


K M PP+K p p + K p 


(137) 


Pitch 


l yA + (4- - 4 ) rp + 4 ( qp - r ) + m [z CG (u - vr + wq ) - x CG ( w - uq + vp)] = M HS + 

M w\ w \ w \ w \ + M q \ q f} kl++ M q q + m ,/p + M m m + M uw uw + M vp v P + M ms p u 15 P ~ ( : 13 8 ) 
M 


6. Yaw 

4'" + (4 - 4)M + 4 (? + r 4 + w [*cg (v- W P + wr) - v CG (m - vr + wq)] = N j 


HS 


+ 


(139) 


(140) 


4|v| V M + N M r H + N *V + + N “' MV + N P« Pq + N ‘" Ur + N «P WP + N «*r U Sr - N P ro P 

Combining acceleration terms invokes an element of convenience for further 
manipulation into matrix form: 

(m -X .) u + mz CG q - my CG r = X HS + X^u \u\ + [X wq - m)wq + (X qq - mx CG )q 2 

+ ( X w + m ) vr + ( x rr + mx co )r 2 ~ my CG pq - mz CG pr + X 
(m -Y v )v-mz CG p + (mx CG - Y .)f = Y m + 7 |v| v|v| + 7 |r| r |r| + my CG r 2 + (Y r -m)ur 

+ (4 - m ) W P + ( Y Pq ~ mx CG ) pq + Y uv uv + my CG p 2 (1 ■ 41 ) 
+mz CG qr + Y mS u 1 8 r + Y prop 

(m ~ Z.) w+my CG p ~{mx CG +Z.)q = Z HS + Z^w|w| + Z^q \q\ + (Z uq + m)uq 

+ (Z, p - m) vp + (4 - mx C G) >P + z m uw (142) 

+mz CG (p 2 +q 2 )~my CG rq + Z uuS u 2 8 p + Z 

~ mz ccX + m .y CG w +(4 ~ k p)p= k hs + K P \ P p\ p\ - (4 - 4 ) qr ~ 4 ( r2 ~ ql ) 

+my CG (; uq-vp)-mz CG ( wp-ur) + K prop 


prop 


(143) 
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mz Cg u - (mx CG + M)w + [l yy -M.)q + I yz r = M HS + M^w\w\ + M^q\q\ 

+ ( M uq -mx CG )uq + (m vp - mx CG ) vp 

+ [ M rp ~ (4 - 4 )] W - lyzQP ( 144 > 

+mz CG (vr -wq) + M uw uw 

+M mS u 2 S p+ M prop 


)ur 


-Wcg «+ { mx cG - 4) v+ iyA +(4 - 4)= 4® + 4|v| v H+4|4 H+(4- - mx CG ) 1 

+ {N wp +mx CG )wp + [N pq -(4 -4)]m (145) 
+I yz rp - my CG (vr - wq ) + N m ,MV + N uuS u 2 8 r 

+N pwp 

Combining the coefficients and variables on the right-hand side of each equation into a 
j^y,/ 7 vector results in the following matrix form: 


(m~ 4 ) 0 

0 (m-n) 

0 0 
0 -fflz 


fflZ, 


CG 


CG 


0 


0 

0 

my CG 

(- mx CG -M w ) 


~my CG (mx CG -N v ) 


0 

mz CG 

-my CG 




~mz CG 

0 

(mx CG -Y r ) 

U 

V 


Z* 

Wcg 

(-mx CG -Z.) 

0 

w 


Z* 


0 

0 

p 


Z* 

0 

(I -M.) 

-I 

q 


Y. M 

0 

\ yy q ) 

-4 

yz 

(4-V) _ 

r 




(146) 


where the mass matrix is invertible, enabling solution of the acceleration vector directly 
using MATLAB/SIMULINK software. 


C. IMPLEMENTATION INTO SIMULINK 

The preceding matrix equations were implemented into the following SIMULINK 

system: 
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Figure 17. SIMULINK Block Diagram of the REMUS Vehicle 


This vehicle model requires the following two inputs: 

• Vehicle Initial Conditions - the initial state vector: 
x = [u V w p q r Xpos y v 


pos y pos pos 


<t> 0 v\ 


where u, v, w are measured in 


m/ 


p, q, r are measured in 


deg/' 


(147) 
x, y, z are 


measured in [ml, and (/), 6 , (//arc measured in [deg]. 


• Vehicle Control Inputs - including rudder and plane pitch angles and thruster and 
main propeller rotation speeds: 


u 


=[4 < 5 , 


p n FLT n SLT n FVT n SlT W 


(148) 
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where S r , 8 p are rudder and plane commands in [deg], n FLT , n SIT , n FVT , n SVT , n are all 
thruster inputs in [RPM ]. 

The input vectors are executed in a feedback loop through the Embedded 
MATLAB Function remus, which contains all vehicle parameters and coefficients and 
solves the following equation: 


. —1 


1 

1 _ 

0 

0 

0 

mz CG 

-my CG 

u 

V 


0 


0 

-mz CG 

0 

(mx CG -Y 

IV 


0 

0 

0 «- z „) 

Wcg 

(-mx CG -Z.) 

0 

p 


0 

-mz CG 

my CG 


0 

0 

q 


mz CG 

0 

{-nvc CG -M v ) 

0 

( 4 -^) 

-4 

r 







~Wcg 
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I* 

2 * 

2 * 
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where the acceleration vector is integrated directly, thereby obtaining the output state 


vector, as well as the output body force and thruster force vectors and total vehicle speed. 
The code for the Embedded MATLAB Function remus is presented in its entirety in 
Appendix C. 


D. VEHICLE SIMULATIONS 

The REMUS vehicle model was simulated using the following control vector 

input: 

u = [ 0 0 2000 -2000 0 0 of (150) 

resulting in the following vehicle motions: 
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Z Direction [m] 



-Remus Vehicle Path 


Remus Vehicle 3D Position Plot 


X Direction [m] 


Figure 18. 3-D Position Plot, Standard Mission 


Remus Vehicle XY Plane Position Plot 



Figure 19. X-Y Plane Position Plot, Standard Mission 
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Remus Vehicle XZ Plane Position Plot 



Figure 20. X-Z Plane Position Plot, Standard Mission 


Remus Vehicle YZ Plane Position Plot 



Figure 21. Y-Z Plane Position Plot, Standard Mission 
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Clearly, the vehicle is inherently unstable. The location of the cross-body 
thrusters relative to the center-point of the vehicle causes an unequal series of forces and 
moments about the body if equal thrust commands are applied, resulting in the sinusoidal 
oscillations along the path of advance. While not directly input, the path of advance is a 
result of these unequal force and moment distributions. Furthermore, in accordance with 
the standard North-East-Down frame of reference, the vehicle has a tendency to rise to 
the surface, denoted in the increasing negative values along the z-axis. This is a direct 
result of the inherent positive buoyancy of the REMUS vehicle in a standard operating 
environment. 

Forcing the vehicle into neutral buoyancy, W - B , results in the following vehicle 
position plots: 


Remus Vehicle 3D Position Plot 



Figure 22. 3-D Position Plot, Neutral Buoyancy 
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Remus Vehicle XY Plane Position Plot 



Figure 23. X-Y Plane, Neutral Buoyancy 


Remus Vehicle XZ Plane Position Plot 



Figure 24. X-Z Plane, Neutral Buoyancy 
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Remus Vehicle TZ Plane Position Plot 



Figure 25. Y-Z Plane, Neutral Buoyancy 

The oscillatory behavior remains; however, the vehicle tends to dive. This tendency is 
attributed to vehicle instability in roll. Implementation of a series of control algorithms is 
required to provide complete vehicle stability. 
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VI. CONTROL IMPLEMENTATION 


A. INTRODUCTION 

This chapter presents the implementation of a Proportional-Integral-Derivative 
(PID) controller to stabilize motions in yaw rate, r. PID controllers are widely used 
because they are generally applicable to most control systems and are particularly useful 
when a complete mathematical model of the plant is unknown [31]. 

The state space and thruster control models developed in Chapters IV and V 
provide a first order differentiable system resulting in simulation estimates on required 
lateral thruster inputs to minimize vehicle motions in the x-y plane. These estimates are 
then incorporated into the PID controller. 

This chapter then discusses the validity of control all motions in the x-y plane, 
u, v, r, i//, x pos , y pos , utilizing only the lateral cross-body thrusters to force rotation about 

the z-axis without deviation in the x-y plane. 

B. PID CONTROL OF YAW RATE 

A PID controller is a feedback mechanism that calculates the error between a 
desired operating point and the measured value. For yaw rate, if/ , the error equation is: 

(151) 

where ij/ e (?) is the error in yaw rate as a function of time, yf c is the constant commanded 

yaw rate and ys(f) is the measured yaw rate as a function of time from the REMUS 

computer simulation. The general block diagram for a PID controller is presented in 
Figure 26. 
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Student Version of MATLAB 


Figure 26. Yaw Rate PID Controller Block Diagram 


In Figure 26, K p , K n K D are the proportional, integral and derivative gains, 
respectively. The values shown in Figure 26 were chosen via trial and error. In fact, 
what is shown in the figure is a Proportional-Derivative (PD) controller rather than a full 
PID controller. The input to the P, I, and D blocks is the difference between the 
commanded angular velocity and the angular velocity calculated through the REMUS 
vehicle model (i.e. the output of Equation (151)). The commanded angular velocity was 


set to y/ c = 7.0 


deg/' 


, with the identical control input vector from Chapter V Section D 


and standard mission profile of slight positive buoyancy. Figures 27 through 31 present 
the results of simulation with a closed yaw-rate loop: 
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Z Direction [m] 


Remus Vehicle 3D Position Plot 



-Remus Vehicle Path 


X Direction [m] 


Figure 27. 3-D Position Plot, PID Control Implemented 


Remus Vehicle XY Plane Position Plot 



Figure 28. X-Y Plane Position Plot, PID Control Implemented 
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Remus Vehicle XZ Plane Position Plot 



Figure 29. X-Z Plane Position Plot, PID Control Implemented 



Figure 30. Pleading Plot, PID Control Implemented 
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Yaw Rate PID Control 



Figure 31. PID Control for Yaw Rate 


As expected, the PID controller forces the yaw rate to converge to 


commanded rate of 7 


deg/' 


, resulting in the linear heading change shown in Figure 


the 

30. 


The positive buoyancy still results in the vehicle gradually surfacing, but the severe 
oscillatory behavior from inherent vehicle instabilities subsides. 


C. MOTION CONTROL IN THE X- Y PLANE 

Restricting vehicle control input to the forward and aft lateral cross-body thrusters 
only, motion control using a linearized model of the REMUS vehicle model in the 
horizontal plane is examined. 


1. Model Linearization - Coupled Thruster Control 

Using the following state vector: 

x = [u v r i//] T (152) 


assuming the forward and aft lateral thrusters are coupled: 
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A system is completely controllable when the controllability matrix possesses full 
rank, where the controllability matrix is defined as: 

Co = [b AB A 2 B 4 3 b] (156) 

Substituting the appropriate coefficient values found in Appendix B yields the following 
matrices: 

-0.0468 0 -1.0490 Ol [ 0 " 

0 -1.9807 -1.3368 0 0.0399 

,B = 

0 0.0307 -35.2346 0 0.0016 

0 0 10 0 

(157) 

0 -0.0017 0.0582 -2.0467 

0.0399 -0.0812 0.2350 -3.0702 

C = 

0.0016 -0.0554 1.9486 -68.6513 

0 0.0016 -0.0554 1.9486 

where n was input as 2,000 RPM. The rank of the controllability matrix is four, and 

therefore full, indicating the system is controllable. Implementation of this controllable 
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linear model indicates that the vehicle can be controlled via commanded states, such as a 
desired heading, yJ command , or desired location in the x-y plane. 

2. Model Linearization - Thruster Differential Analysis 

Validation of the non-linear model requires an analysis of the applied moment in 
both the horizontal and vertical planes of motion. This differential is then implemented 
into the REMUS vehicle itself for comparison to the computer model. 

a. Lateral Thruster Differential 

Due to non-symmetrical placement of the cross-body thrusters, an inherent 
moment exists when equal thruster RPM commands are input into the model. Using the 
linearized model, several simulations were run holding the forward thruster constant and 
varying stern lateral thruster speeds analyzing vehicle response to a commanded yaw rate 

of 5 . Table 10 presents the results: 


FLT [RPM] 

SLT [RPM] 

Thruster 

Differential [RPM] 

500 

-300 

200 

1000 

-700 

300 

1500 

-1400 

100 

2000 

-1900 

100 

2500 

-2300 

200 

3000 

-2600 

400 

3500 

-3300 

200 

4000 

-3600 

400 

4500 

-4000 

500 

5000 

-4500 

500 


Table 10. Lateral Cross-Body Thruster Differential 
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This results in an average thruster differential of 290 RPM. Figure 32 presents the 
implementation of this differential into the linear model with no PID control: 


Remus Vehicle XY Plane Position Plot 



n . - w n i • -1 - n h i •_» ji - -i - i 






X Direction [m] 

Figure 32. X-Y Plane, Linear Model 

Applying the constant thruster differential to the model minimizes vehicle motions in the 
x-y plane. 


b. Vertical Thruster Differential 

The vertical thruster differential was determined using a similar analysis 
applied to the lateral thrusters. Once more, the forward vertical thruster was held 
constant and the stem vertical thruster was varied to achieve a commanded pitch of 
0[deg]. Table 11 presents the results: 
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FVT [RPM] 

SVT [RPM] 

Thruster 

Differential [RPM] 

500 

-200 

300 

1000 

-550 

450 

1500 

-1100 

400 

2000 

-1550 

450 

2500 

-2100 

400 

3000 

-2550 

450 

3500 

-2950 

550 

4000 

-3400 

600 

4500 

-4000 

500 

5000 

-4450 

550 


Table 11. Vertical Cross-Body Thruster Differential 


These simulations result in an average vertical thruster differential of 465 RPM. Figure 
33 presents the implementation of the averaged thruster differential into the linear model 
with no PID control. 


Remus Vehicle Pitch Plot 



Figure 33. Vehicle Pitch, Linear Model 
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The thrust differential briefly forces the vehicle to small perturbations in pitch, 
but requires detailed controllability analysis and implementation of a controller to prevent 
the instability in the vehicle model from growing unbounded as in Figure 33. 


68 



VII. CONCLUSIONS 


A. INTRODUCTION 

This chapter presents an analysis of results, conclusions drawn and possible 
avenues of future study. Recall that non-linear vehicle model simulations utilized 
equivalent thruster control to force the vehicle to rotate about its center, minimizing 
vehicle movement in the horizontal plane, while the linear vehicle model analyzed the 
thruster induced moments on the vehicle. 

B. ANALYSIS OF RESULTS 

1. Hydrodynamic Coefficients 

A significant amount of time was spent developing the required hydrodynamic 
coefficients for this latest version of the REMUS vehicle. These coefficients As stated 
previously, the coefficients remain dimensional for ease of comparison with the 
coefficients from the original REMUS vehicle [12]. Significant deviations in coefficient 
values occurred. A simple unit analysis for the non-linear force coefficients shows that 
the primary contributions to many of the coefficients are from mass [ kg | and a 

length [m], either radial or axial. Given that the vehicle mass and axial length is 
significantly larger; their impact on respective coefficients is proportionately larger. 
Similarly, the non-linear moment coefficients generally possess units of mass [kg\ and 

area [ m 2 ], resulting in even larger deviations from values detennined in [12]. 

2. Vehicle Modeling 

Implementation of the derived hydrodynamic coefficients results in an accurate 

non-linear model of true REMUS vehicle motions. This is evident in the required 

thruster differential due to non-symmetrical placement of the cross-body thrusters. 

Simulations using the non-linear model with equal and opposite thruster inputs results in 

a sinusoidal path of advance caused by unequal moments generated by the thruster pairs. 

To minimize these motions in the horizontal plane, the stem lateral thruster requires a 

setting of 290 RPM less than the forward lateral thruster to eliminate the unequal 
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moment. Similarly, for the vertical thrusters, a differential of 465 RPM is required to 
minimize vehicle pitch while ascending or descending. However, a more thorough 
comparison between actual in-water data and simulation results needs to be conducted to 
fully validate this REMUS vehicle model. Tow tank experiments validating the radial 
and axial drag coefficients are of particular interest, since the significant size differential 
forces large differences between the derived hydrodynamic terms in this work and those 
in [12]. 


3. Vehicle Control 


Implementation of a P1D control for yaw rate through controllability analysis of a 
linear model further supports the validity of the REMUS vehicle model. Although only a 
single batch of simulation results is presented; the vehicle model tracks as commanded to 


any desired yaw rate within the model boundaries of 0 


deg. 


to 14.5 


deg. 


. The yaw 


rate boundary limits are dependent upon the minimum and maximum thruster range of 0 
RPM to 5,000 RPM. In depth thruster analysis remains to be executed to provide a fully 
controllable thruster model in closed loop control. Of particular interest is the 
contribution of individual thruster pairs to their respective angular velocities and 
accelerations, as well as the effects of slaved and non-slaved thruster pairs to the vehicle 
dynamics. Implementation of closed loop control requires this detailed analysis of the 
thrusters on the vehicle dynamics and would facilitate faster response time to track the 
desired yaw rate. 


C. FUTURE WORK 

The modularity of the REMUS vehicle provides a unique foundation for 
continuous work in the field of autonomous underwater vehicles. This work provides a 
logical origin for follow-on study, particularly in the following fields: 

• Object Avoidance - Implementation of cross-body thruster control during object 
avoidance maneuvers can be analyzed and compared to results found in [13], 
enabling further study in the field of optimal control, generating logic algorithms 
optimizing thruster and control surface use to achieve a goal or series of goals. 
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• Simultaneous Localization and Mapping (SLAM) - The robust model reduces 
errors in state space estimates, thereby reducing the area of uncertainty in the 
localization process. 

• Model Reference Adaptive Control (MRAC) - A valid vehicle model enables the 
implementation of MRAC to provide adjustable closed loop control of the 
REMUS vehicle. 

• Multiple AUV Operations - A robust model of each vehicle in simultaneous 
operations allows more coordinated command and control by the designated lead 
vehicle to ensure reliable operation and mission completion. 

• AUV Operations in Larger Operating Envelopes - A robust on-board model 
permits the AUV to operate outside the normal environmental tolerances. 
Feedback sensors combined with sea-state prediction algorithms would expand 
AUV use to encompass a wider array of military and civilian operating 
environments. 
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APPENDIX A. TABLES OF PARAMETERS 


A.l STANDARD REMUS HULL PARAMETERS 


Parameter 

Value 

Units 

Description 

P 

1.03e+03 

P 

1-1 

Seawater Density 

L 

2.72415e+00 


H 


Vehicle Length 

D 

1.9e-01 


H 


Vehicle Maximum Diameter 

m 

6.60399e+l 


U 


Vehicle Mass 

W 

6.4785e+2 




Vehicle Weight 

B 

6.9799e+2 


[N] 


Vehicle Buoyancy 

A f 

2.84e-2 

[ 

m 2 

] 

Hull frontal area 

A P 

4.786e-l 

[ 

m 2 

} 

Hull projected area (x-z plane) 


Table 12. REMUS Hull Parameters 


A.2 CENTER OF BUOYANCY RELATIVE TO ORIGIN AT VEHICLE NOSE 


Parameter 

Value 

Units 

Description 

%CB 

-1.362075e0 

H 

Longitudinal Center of Buoyancy 

VCB 

0.00 

H 

Transverse Center of Buoyancy 

ZCB 

0.00 

H 

Vertical Center of Buoyancy 


Table 13. Center of Buoyancy 
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A.3 CENTER OF GRAVITY RELATIVE TO ORIGIN AT VEHICLE HALF 
LENGTH 


Parameter 

Value 

Units 

Description 

*CG 

1.4325e-2 

H 

Longitudinal Center of Gravity 

ycc 

0.00 

H 

Transverse Center of Gravity 

ZCG 

-4.7625e-2 

H 

Vertical Center of Gravity 


Table 14. Center of Gravity 

A.4 HULL COORDINATES FOR LIMITS OF INTEGRATION RELATIVE 
TO ORIGIN AT VEHICLE HALF LENGTH 


Following the notation used in [12]. 


Parameter 

Value 

Units 

Description 

A 

-1.362075e0 

H 

Aft End of Tail Section 


-6.40975e-l 

H 

Forward End of Tail Section 

A/ 

-1.235075e0 

H 

Aft End of Fin Section 

X h 

-1.150975e0 

H 

Forward End of Fin section 

A 

8.35025e-l 

H 

Aft End of Bow Section 


1.360275e0 

H 

Forward End of Bow Section 


Table 15. Hull Coordinates for Limits of Integration 
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A.5 STANDARD REMUS FIN PARAMETERS 


Parameter 

Value 

Units 

Description 

S fin 

6.65e-03 

[ra ! ] 

Plan fonn area 

bfin 

8.57e-02 

H 

Fin Span 

a fin 

1.39e-01 

H 

Maximum Fin Height above 
centerline 


1.17e-01 

H 

Average Fin Height above 
centerline 

C 4 

5.58e-01 


Cross-flow drag coefficient 

t 

5.31 e-01 


Fin Taper Ratio 

X fin 

-6.85e-01 

H 

Fin post location relative to origin 
/ 

at — 

2 


Table 16. Standard REMUS Fin Parameters 
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APPENDIX B. TABLES OF COMBINED NON-LINEAR 

COEFFICIENTS 


B.l NON-LINEAR FORCE COEFFICIENTS 


Parameter 

Value 

Units 

Description 


-3.13e0 


kg/ 

/ m 


Axial Drag coefficient in surge 


-8.84e-l 

[**] 

Axial Added Mass coefficient in 
surge 

K 

-7.78e+l 


kg/ 

/ rad 


Surge cross-term in heave and 
pitch 

V 

-4.16e0 


kg ■ m / 

/ rad 


Surge drag coefficient in pitch 

x vr 

7.78e+l 

■ 

kg/ 

/ rad 


Surge cross-term in sway and yaw 

x rr 

-4.16e0 


kg ■ m / 

/rad 


Surge drag coefficient in yaw 

y w 

-2.85e+2 


kg/ 

/ m 


Axial drag coefficient in sway 

Y rr 

1.30e+l 


kg ■ m / 

/ rad 2 _ 


Sway coefficient in yaw 

K v 

1.16e+l 


kg/ 

/ m 


Body lift force and fin lift 

n 

-7.78e+l 

[**] 

Axial Added Mass coefficient in 
sway 

Y> 

4.16e0 


kg ■ m / 

/ rad 


Added mass coefficient in yaw 

Y ur 

-8.84e-l 

■ 

kg/ 

/ rad 


Sway cross-term in surge and yaw 

Y wp 

7.78e+l 

■ 

kg/ 

/ rad 


Sway cross-term in heave and roll 

Y ,„ 

4.16e0 


kg ■ m / 

/ rad 


Sway cross-tenn in roll and pitch 

Y uu s r 

2.13e+l 


kg ■ m / 

/ rad 


Fin lift force coefficient (rudder) 
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-2.85e+2 


kg/ 

/ m 


Axial drag coefficient in heave 

z 

m 

1.30e+l 


kg ■ m / 

/rad 2 _ 


Heave coefficient in pitch 


1.16e+l 


kg/ 

/ m 


Body lift force and tin lift 


-7.78e+l 

[ks] 

Axial Added Mass coefficient in 
heave 


-4.16e+0 


kg ■ m / 

/rad 


Added mass in pitch 

^uq 

8.84e-l 

['%a, 

■ 

Heave cross-term in surge and 
pitch 

A 

-7.78e+l 

\% P 

- 

Heave cross-term in sway and roll 

Z rp 

4.16e0 

[Vra d , 

- 

Heave cross-term in roll and yaw 

uuS s 

-2.13e+l 


1 1 

1 

1_1 


Fin lift force coefficient (stem- 
plane) 


Table 17. Non-Linear Force Coefficients 

B.2 NON-LINEAR MOMENT COEFFICIENTS 


Parameter 

Value 

Units 

Description 

K 

-1.30e-l 


kg ■ m’/ 

/rad 2 


Axial rolling drag 

K P 

-1.04e-2 


kg ■ in 2 / 

/rad 2 


Added mass moment in roll 

M ww 

-8.82e0 


Cross-flow drag in heave 

M 

qq 

-2.81e+3 


kg ■ m 2 / 

/ rad 2 


Cross-flow drag in pitch 

M uw 

-5.48e-l 

[kg] 

Body lift and tin lift 


-4.16e0 

[kg ■ m ] 

Added mass moment in heave 


-4.52e+l 


1-1 

N 

s 

Co 


Added mass moment in pitch 

M uq 

4.16e0 


kg ■ m / 

/rad 


Added mass moment cross-term 
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in surge and pitch 

K 

-4.16e0 


kg ■ m / 

/rad 


Added mass moment cross-term 
in sway and roll 

K 

4.52e+l 

- 

kg ■ m 2 / 

/rad" 

~ 

Added mass moment cross-term 
in roll and pitch 

M uu s r 

-1.46e+l 


Fin lift moment (rudder) 

N vv 

8.82e0 

[**] 

Cross-flow drag in sway 

N rr 

-2.81e+3 

- 

kg ■ m 2 / 

/rad" 

~ 

Cross-flow drag in yaw 

N uv 

5.48e-l 

[**] 

Body lift and fin lift 


4.16e0 

[kg ■ m ] 

Added mass moment in sway 

N r 

-4.52e+l 


1-1 

"73 

N 

s 

Co 


Axial added mass moment in yaw 

N ur 

4.16e0 


kg ■ m / 

/rad 


Added mass moment cross-term 
in surge and yaw 

K P 

-4.16e0 


kg ■ m / 

/rad 


Added mass moment cross-term 
in heave and roll 

N „ 

-4.52e+l 

- 

kg ■ m 2 / 

/rad" 

~ 

Added mass moment cross-term 
in roll and pitch 

^uuS s 

-1.46e+l 

[Va d \ 

Fin lift moment (stern-plane) 


Table 18. Non-Linear Moment Coefficients 
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APPENDIX C. MATLAB CODE 


C.l EMBEDDED MATLAB FUNCTION: REMUS.M 

function [x_dot,p_dot, adot, F, F_prop, U] = remus(x,ui) 

%% Author: Sean Doherty 
%% Date: 30 September 2010 
%% Revision: 17th 
%% Revision Date: 27 February 2011 
%% References : 

% L. R. Fodrea, "Obstacle Avoidance Control for the REMUS Autonomous 
% Underwater Vehicle." M.S. Thesis, Naval Postgraduate School, 

% Monterey, California, USA, 2002. [Print]. 

% 

% T. I. Fossen and T. Perez, Marine Systems Simulator (MSS),_ 

% developed at the Norwegian University of Science and 
% Technology, 2010. [Online] Software Available: 

% http://www.marinecontrol.org/. [Accessed: 5 March 2010]. 

% 

% T. Prestero, "Verification of a Six-Degree of Freedom Simulation 
% Model for the REMUS Underwater Vehicle." M.S. Thesis, 

% Massachusetts Institute of Technology, Cambridge, 

% Massachusetts, USA, 2001. [Print]. 

% 

%% Details: 

% 

% [xdot,F,U, F_prop] = remus(x,ui) gives the derivative of the state vector: 

%x = [uvwpqrxyz phi theta psi ]' for the REMUS AUV. All 
% requisite parameters and hydrodynamic coefficients are included within 
% this function. 

% _ 

% Inputs: 

% The state vector is defined as: 

% 

% x = [u v w p q r xpos ypos zpos phi theta psi]' where 

% 

% u = surge velocity (m/s) 

% v = sway velocity (m/s) 

% w = heave velocity (m/s) 

% p = roll velocity (rad/s) 

% q = pitch velocity (rad/s) 

% r = yaw velocity (rad/s) 

% xpos = position in x-direction (m) 

% ypos = position in y-direction (m) 


81 



% zpos = position in z-direction (m) 

% phi = roll angle (rad) 

% theta = pitch angle (rad) 

% psi = yaw angle (rad) 

% 

% The input vector is : 

% 

% ui = [ delta_r delta_p fit sit fvt svt n ]' where 

% 

% delta_r = rudder angle (deg) 

% delta_p = port and starboard stern plane (deg) 

% fit = forward lateral thruster (RPM) 

% sit = stern lateral thruster (RPM) 

% fvt = forward vertical thruster (RPM) 

% svt = stem vertical thruster (RPM) 

% n = main propeller shaft speed (percent max thrust) 

% 

% _ 

% Outputs: 

% 

% Speed: U = (m/s) 

% Forces: F =[XYZKMN]' 

% F_prop = [X_prop Y_prop Z_prop K_prop M_prop N_prop]’ 

% Accelerations: xdot = [Xdot Ydot Zdot Kdot Mdot Ndot]' 

% 

% _ 

%% Note: 

% Portions of this code were adapted from the MSS GNC toolkit. 

% The MSS GNC is a Matlab toolbox for guidance, navigation and control. 

% The toolbox is part of the Marine Systems Simulator (MSS). 

% 

% Copyright (C) 2008 Thor I. Fossen and Tristan Perez 

% 

% This program is free software: you can redistribute it and/or modify 
% it under the terms of the GNU General Public License as published by 
% the Free Software Foundation, either version 3 of the License, or 
% (at your option) any later version. 

% 

% This program is distributed in the hope that it will be useful, but 
% WITHOUT ANY WARRANTY; without even the implied warranty of 
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
% GNU General Public License for more details. 

% 

% You should have received a copy of the GNU General Public License 
% along with this program. If not, see <http://www.gnu.org/licenses/>. 
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% 

% E-mail: contact@marinecontrol.org 
% URL: <http://www.marinecontrol.org> 

% _ 

%% Verification of input and state dimensions 

if (length(x) ~= 12),error('x-vector must have dimension 12 !');end 

if (length(ui) ~= 7),error('u-vector must have dimension 7 !');end 

%% Dimensional states 
u = x(l); v =x(2); w = x(3); 
p =x(4); q = x(5); r = x(6); 
phi = x(10); theta = x(l 1); psi = x(12); 

U = sqrt(u A 2+v A 2+w A 2); % speed 

% Reorganizing dimensional states: 

%nu = [u v w p q r]’; 

%% REMUS Vehicle Measurements (taken 28 October 2010): 

%L = 2.72415; % Length in meters (m) [includes BlueView Package] 

%D = 0.19; % Maximum diameter in meters (m) 

g = 9.81; % Acceleration due to gravity in meters/second A 2 (m/s A 2) 

m = 66.0399; % Mass in kilgrams (kg) 

W = m*g; % Weight in Newtons (N) 

B = (m+0.5)*g; % Measured Vehicle Buoyancy (N) 

V_max = 2.881; % Maximum Forward Velocity (m/s) 

% Moments of Inertia WRT Origin at Half-Length 

I_xx = 1.154E0; % kg*m A 2 

I_yy = 3.470E1; % kg*m A 2 

I_zz = I_yy; % kg*m A 2 

I_yz = 5.768E-1; 

I_zy = I_yz; 


% Center of Buoyancy WRT Origin at Vehicle Nose 
x_cb = -1.36208E0; % x-location (m) 

y_cb = 0.00; % y-location (m) 

z_cb = 0.00; % z-location (m) 

% Center of Gravity WRT Origin at Vehicle Nose 
x_cg = 1.4325E-2; % x-location (m) 

y_cg = 0.00; % y-location (m) 

z_cg =-1.016E-1; % z-location (m) 
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% Non-Linear Force Coefficients 


X_uu = -3.13E0; 
X_udot = -8.84E-1; 
X_wq = -7.78E1; 
X_qq = -4.16E0; 
X_vr = 7.78E1; 
X_rr = -4.16E0; 


% Cross-flow Drag (kg/m) 

% Added Mass (kg) 

% Added Mass Cross-term (kg/rad) 

% Added Mass Cross-term (kg*m/rad) 
% Added Mass Cross-term (kg/rad) 

% Added Mass Cross-term (kg*m/rad) 


Y_vv = -2.85E2; 
Y_rr = -1.30E1; 
Y_uv = 1.16E1; 
Y_vdot = -X_vr; 
Yrdot = -X_rr; 
Y_ur = Xudot; 
Y_wp = -X_wq; 
Y_pq = -X_qq; 
Yuudr = 2.13E1; 


% Cross-flow Drag (kg/m) 

% Cross-flow Drag (kg*m/rad A 2) 

% Body Lift Force and Fin Lift (kg/m) 

% Added Mass (kg) 

% Added Mass (kg*m/rad) 

% Added Mass Cross-term and Fin Lift (kg/rad) 
% Added Mass Cross-temi (kg/rad) 

% Added Mass Cross-term (kg*m/rad) 

% Fin Lift Force (kg/(m*rad)) 


Z_ww = Y_vv; 
Z_qq = -Y_rr; 
Z_uw = Y_uv; 
Z_wdot= X_wq; 
Zqdot = X_qq; 
Z_uq = -Xudot; 
Z_vp = Y_vdot; 
Z_rp = Yrdot; 
Z_uuds = -2.13E1; 


% Cross-flow Drag (kg/m) 

% Cross-flow Drag (kg*m/rad) 

% Body Lift Force and Fin Lift (kg/m) 

% Added Mass (kg) 

% Added Mass (kg*m/rad) 

% Added Mass Cross-term and Fin Lift (kg/rad) 
% Added Mass Cross-term (kg/rad) 

% Added Mass Cross-term (kg/rad) 

% Fin Lift Force (kg/(m*rad)) 


% Non-Linear Moment Coefficients 


K_pp = -1.30E-1; % Rolling Resistance (kg*m ' 2/rad / 2) 

K_pdot = -1.04E-2; % Added Mass (kg*m A 2/rad) 


M_ww = -8.82E0; 
M_qq = -2.81E3; 
M_uw = -5.48E-1; 
Mwdot = Zqdot; 
M_qdot = -4.52E1; 
M_uq = -X_qq; 
M_vp = -M_uq; 
M_rp = 4.52E1; 
Muuds = -1.46E1; 


% Cross-flow Drag (kg) 

% Cross-flow Drag (kg*m A 2/rad A 2) 

% Body and Fin Lift and Munk Moment (kg) 

% Added Mass (kg*m) 

% Added Mass (kg*m A 2/rad) 

% Added Mass Cross-term and Fin Lift (kg*m/rad) 
% Added Mass Cross-term (kg*m/rad) 

% Added Mass Cross-term (kg*m A 2/rad A 2) 

% Fin Lift Moment (kg/rad) 


N_vv = -M_ww; % Cross-flow Drag (kg) 
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N_rr = M_qq; 
N_uv = 5.48E-1; 
N_vdot= Y_rdot; 
N_rdot = M_qdot; 
N_ur = M_uq; 
N_wp = -N_ur; 
N_pq = -M_rp; 
N_uudr= M_uuds; 


% Cross-flow Drag (kg*m A 2/rad A 2) 

% Body and Fin Lift and Munk Moment (kg) 

% Added Mass (kg*m) 

% Added Mass (kg*m A 2/rad) 

% Added Mass Cross-term and Fin Lift (kg*m/rad) 
% Added Mass Cross-term (kg*m/rad) 

% Added Mass Cross-term (kg*m A 2/rad A 2) 

% Fin Lift Moment (kg/rad) 


% Rigid Body Mass Matrix: 


Mil = [m0 0; 
0 m 0; 

0 0 m]; 


M12 = [ 0 (m*z_cg) (-m*y_cg); 

(-m*z_cg) 0 (m*x_cg); 

(m*y_cg) (-m*x_cg) 0 ]; 

M21 = [ 0 (-m*z_cg) (m*y_cg); 

(m*z_cg) 0 (-m*x_cg); 

(-m*y_cg) (m*x_cg) 0 ]; 


M22 = [I_xx 0 0; 

o Lyy -i_yz; 

0 -I_zy I_zz]; 

% Added Mass Matrix: 
MAI 1 = [X_udot 0 0; 

0 Y_vdot 0; 

0 0 Z wdot]; 


MA12 = [ 0 0 0; 

0 0 Y rdot; 

0 Z qdot 0]; 


MA21 = [ 0 0 0; 

0 0 M_wdot; 

0 N_vdot 0]; 

MA22 = [K_pdot 0 0; 

0 M qdot 0; 

0 0 N_rdot]; 

% Total Mass Matrix (M-MA): 
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MM = [(Ml 1-MA11) (M12-MA12); 
(M21-MA21) (M22-MA22)]; 

InvM = inv(MM); 


%% 


cl = cos(phi); 
c2 = cos(theta); 
c3 = cos(psi); 
si = sin(phi); 
s2 = sin(theta); 
s3 = sin(psi); 
t2 = tan( theta); 

%% Control input (rudder, plane, thrusters and propeller) 
delta r = ui(l); % Angle (rad) 

delta_p = ui(2); % Angle (rad) 

% Thruster Terms (Cross-body thruster data taken 1 November 2010) 

% Main Thruster (from Prestero Thesis): 

n_max = X_uu*V_max*abs(V_max); % Propeller Max Thrust (N) 

n_maxT = (0.995*((W*y_cg)-(B*y_cb)))-(0.093*((W*z_cg)-(B*z_cb))); % Propeller 

Max Torque (N*m) 

n = ((ui(7)/100)*n_max); % Thrust (N) 

n_t = ((ui(7)/100)*n_maxT); % Torque (N-m) 

% Cross-body Thruster Thrust (N): 

fit = (-4E-14*ui(3) A 4)+(3E-10*ui(3) A 3)-(2E-7*ui(3) A 2)+... 

(0.0004*ui(3))-0.0054; 

sit = (-4E-14*ui(4) A 4)+(3E-10*ui(4) A 3)-(2E-7*ui(4) A 2)+... 

(0.0004*ui(4))-0.0054; 

fvt = (-4E-14*ui(5) A 4)+(3E-10*ui(5) A 3)-(2E-7*ui(5) A 2)+... 

(0.0004*ui(5))-0.0054; 

svt = (-4E- 14*ui(6) A 4)+(3E-10*ui(6) A 3)-(2E-7*ui(6) A 2)+... 

(0.0004*ui(6))-0.0054; 

% Cross-Body Thruster Torque (N-m): 

fit t = (-4.702E-14*ui(3) A 4)+(3.302E-10*ui(3) A 3)-(2.700E-7*ui(3) A 2)+... 
(4.360E-4*ui(3))-(6.290E-3); 

sit_t = (1.653E-14*ui(4) A 4)-(1.161E-10*ui(4) A 3)+(9.489E-8*ui(4) A 2)-... 
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(1.532E-4*ui(4))+(2.21 IE-3); 

fvt_t = (-4.523E-14*ui(5) A 4)+(3.176E-10*ui(5) A 3)-(2.597E-7*ui(5) A 2)+... 
(4.193E-4*ui(5))-(6.050E-3); 

svt_t = (E832E-14*ui(6) A 4)-(1.287E-10*ui(6) A 3)+(1.052E-7*ui(6) A 2)-... 

(1.699E-4*ui(6))+(2.45 IE-3); 

%% Cross-body Thruster Forces and Moments: 

X_prop = n; 

Y_prop = fit + sit; 

Z_prop = fvt + svt; 

K_prop = n_t; 

M_prop = flt_t + sltt; 

N_prop = fvt_t + svt_t; 

F_prop = [X_prop Y_prop Z_prop K_prop M_prop N_prop]'; 

%% Total Forces and Moments from Equations of Motion 
% Equations of Motion obtained from the Prestero thesis. Specifically, 

% equations (3.8) and (4.49) are used. The initial time derivatives of 
% the state vector are assumed to be zero 

% Surge Equation of Motion: 

X = -((W-B)*sin(theta)) + (X_uu*u*abs(u)) + ((X_wq-m)*w*q) + ... 

((X_qq-(m*x_cg))*q A 2) + ((X_vr+m)*v*r) + ((X_rr+(m*x_cg))*r A 2) -... 
(m*y_cg*p*q) - (m*z_cg*p*r) + X_prop; 

% Sway Equation of Motion: 

Y = ((W-B)*cos(theta)*sin(phi)) + (Y_vv*v*abs(v)) + (Y_rr*r*abs(r)) + ... 
(m*y_cg*r A 2) + ((Y_ur-m)*u*r) + ((Y_wp-m)*w*p) + ... 
((Y_pq-(m*x_cg))*p*q) + (Y_uv*u*v) + (m*y_cg*p A 2) + ... 
((m*z_cg)*q*r) + (Y_uudr*u A 2*delta_r) + Y_prop; 

% Heave Equation of Motion: 

Z = ((W-B)*cos(theta)*cos(phi)) + (Z_ww*w*abs(w)) + (Z_qq*q*abs(q)) + ... 
((Z_uq+m)*u*q) + ((Z_vp-m)*v*p) + ((Z_rp-(m*x_cg))*r*p) + ... 
(Z_uw*u*w) + ((m*z_cg)*(p A 2+q A 2)) - (m*y_cg*r*q) + ... 
(Z_uuds*u A 2*delta_p) + Z_prop; 

% Roll Equation of Motion: 

K = -(((y_cg*W)-(y_cb*B))*cos(theta)*cos(phi)) - ... 

(((z_cg*W)-(z_cb*B))*cos(thcta)*sin(phi)) + (K_pp*p*abs(p)) -... 
((I_zz-I_yy)*q*r) - (I_yz*(r A 2-q A 2)) + (m*y_cg*((u*q)-(v*p))) -... 
((m*z_cg)*((w*p)-(u*r))) + K_prop; 
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% Pitch Equation of Motion: 

M = -(((z_cg*W)-(z_cb*B))*sin(theta)) -... 

(((x_cg*W)-(x_cb*B))*cos(theta)*cos(phi)) + (M_ww*w*abs(w)) + ... 
(M_qq*q*abs(q)) + ((M_uq-(m*x_cg))*u*q) + ((M_vp-(m*x_cg))*v*p) + ... 
((M_rp-(I_xx-I_zz))*r*p) - (I_yz*q*p) + ((m*z_cg)*((v*r)-(w*q))) + ... 

(M_uw*u*w) + (M_uuds*u A 2*delta_p) + M_prop; 

% Yaw Equation of Motion: 

N = -(((x_cg*W)-(x_cb*B))*cos(thcta)*sin(phi)) + ... 

(((y_cg*W)-(y_cb*B))*sin(thcta)) + (N_vv*v*abs(v)) + ... 

(N_rr*r*abs(r)) + ((N_ur-(m*x_cg))*u*r) + ((N_wp+(m*x_cg))*w*p) + ... 
((N_pq-(I_yy-I_zz))*p*q) + (I_yz*r*p) - (m*y_cg*((v*r)-(w*q))) + ... 

(N_uv*u*v) + (N_uudr*u A 2*delta_r) + N_prop; 

F = [X YZKMN]'; 

%% Accelerations: 
xdot = ... 

[InvM( 1, l)*X+InvM( 1,2)*Y+InvM( 1,3)*Z+InvM( 1,4)*K+InvM( 1,5)*M+InvM( 1,6)*N; 
InvM(2,l)*X+InvM(2,2)*Y+InvM(2,3)*Z+InvM(2,4)*K+InvM(2,5)*M+InvM(2,6)*N; 
InvM(3,l)*X+InvM(3,2)*Y+InvM(3,3)*Z+InvM(3,4)*K+InvM(3,5)*M+InvM(3,6)*N; 
InvM(4,l)*X+InvM(4,2)*Y+InvM(4,3)*Z+InvM(4,4)*K+InvM(4,5)*M+InvM(4,6)*N; 
InvM(5,l)*X+InvM(5,2)*Y+InvM(5,3)*Z+InvM(5,4)*K+InvM(5,5)*M+InvM(5,6)*N; 
InvM(6,l)*X+InvM(6,2)*Y+InvM(6,3)*Z+InvM(6,4)*K+InvM(6,5)*M+InvM(6,6)*N]; 

pdot = ... 

[(c3*c2*u) + (((c3*s2*sl)-(s3*cl))*v) + (((s3*sl)+(c3*cl*s2))*w); 

(s3*c2*u) + (((cl*c3)+(sl*s2*s3))*v) + (((cl*s2*s3)-(c3*sl))*w); 

(-s2*u) + (c2*sl*v)+ (cl*c2*w)]; 


a dot 


[(p)+ (sl*t2*q) + 

(cl*q) - 
((sl/c2)*q) + 


(cl*t2*r); 

(sl*r); 

((cl/c2)*r)]; 
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C.2 REMUS VEHICLE MODEL BLOCK DIAGRAM 
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